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Section 1 


INTRODUCTION 


Calspan Corporation has conducted research on warm fog properties 
and fog modification concepts under sponsorship of the National Aeronautics 
and Space Administration since 1963. Results of previous investigations, 
which have included theoretical and experimental studies of warm fog modifi- 
cation concepts, numerical simulations of fog seeding experiments, numerical 
modeling studies of radiation and advection fog, and fog life cycle studies 
have been documented in numerous NASA contractor reports (see e.g., NASA 
reports, CR-72, CR-368, CR-675, CR-1071, CR-1731, SP-212, CR-2078, CR-2079, 
and CR-2201). These investigations have played an important part in our 
understanding of fog and the complex meteorological processes that are involved 
in its formation, persistence, and dissipation. This report describes the 
results of the most recent investigation under Contract No. NAS8- 30776 with 
George C. Marshall Space Flight Center. 

Contract performance was under the technical cognizance of Mr. O.H. 
Vaughan, Jr., Aerospace Environment Division, Space Science Laboratory of 
the NASA Marshall Space Flight Center, Huntsville, Alabama. Mr. John Enders, 
Chief of the Aviation Safety Technology Branch, NASA Office of Aeronautics 
and Space Technology, provided encouragement and the necessary support for 
the accomplishment of this work. 

Field experiments to test the concept of modifying warm fogs by 
seeding with sized hygroscopic materials and observational studies on both 
radiation and advection fog have shown that a thorough understanding of com- 
plex natural processes of fog formation, development, and dissipation is 
required to effectively design and implement fog modification techniques. 
Working toward this end, a one-dimensional numerical model of radiation fog 
was developed and tested against comprehensive field observations of valley 
fog (Pilie et al., 1972). To further extend understanding of the processes 


1 



which determine the life cycles of fogs, a preliminary two-dimensional numeri- 
cal model of advection fog was developed by Calspan (Mack et al . , 1972) under 
joint sponsorship of NASA and the Naval Air Systems Command. While the use- 
fulness and basic capabilities of the advection fog model were demonstrated, 
considerable potential for important research with the model remained. 

In this investigation, the capabilities of advection fog model were 
further developed and delineated. The model was improved by including prognostic 
equations for the horizontal wind and by incorporating the influences of both 
the predicted wind shear and temperature gradient upon the turbulent exchange 
coefficients. To establish the range of validity of the model and to delineate 
important areas for further development and application of the model, fog 
formation and development processes over an ocean surface were simulated 
using the advection fog model, and the results were compared with recent 
Calspan observational studies of sea fog (Mack et al . , 1975). The physical 
and mathematical foundations of the advection fog model and the results of 
numerical experiments with the model are discussed in Section 2. Recommenda- 
tions are presented for future modeling research which builds upon the improve- 
ment and testing of the advection fog model in the present study. 

In addition, a survey was carried out of candidate fog or cloud 
physics experiments for the zero-G environment of a shuttle-type spacecraft. 
Emphasis was placed on defining experiments which are simple in nature and 
relevant to fog modification problems. The survey and recommendations for 
zero-G experiments are presented in Section 3. 
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Section 2 


A NUMERICAL MODEL OF ADVECTION FOG 
2.1 INTRODUCTION 

The principal effort under this program was devoted to improving 
a preliminary two-dimensional numerical model of advection fog developed 
under NASA sponsorship (Mack et al . , 1972) and to establishing the range of 
validity of the improved model. The primary improvements to the model were 
the incorporation of prognostic equations for the horizontal wind and the 
utilization of turbulent exchange coefficients which depend on both the pre- 
dicted wind and potential temperature gradients near the surface. The improved 
model predicts the evolution of potential temperature, water vapor content, 
liquid water content, and the horizontal wind in a vertical plane as determined 
by the processes of vertical turbulent transfer and horizontal advection, as 
well as radiative cooling and drop sedimentation in fog. The model is designed 
to simulate the formation, development, or dissipation of advection fog in 
response to fluxes of heat and moisture to or from the atmosphere to the sur- 
face as driven by horizontal discontinuities in the surface temperature. 

The improved advection fog model was applied to simulate fog forma- 
tion and development over an ocean surface, where the boundary condition at 
the surface leads to the transfer of water vapor and heat to or from the model 
atmosphere in response to advection over surface temperature changes. Results 
obtained from numerical experiments with the advection fog model have been 
compared with results from the Calspan observational studies of fog at sea 
(Mack et al., 1975) to delineate the range of validity of the model and to 
suggest fruitful directions for future model improvements and applications. 

The physical and mathematical foundations of the numerical model 
are presented in Section 2.2. Documentation of the computer program for the 
two-dimensional advection fog model is provided in Appendix A. In Section 2.3, 
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results of numerical experiments on the formation and development of advection 
fog are discussed. Recommendations are presented in Section 2,4 for future 
numerical modeling research which builds upon the advection fog model develop- 
ment and testing in the present study. 

2.2 NUMERICAL MODEL 

2.2.1 Major Assumptions 

The following assumptions are adopted in the numerical modeling study 
of advection fog; 

(a) The model is two-dimensional in the x-z plane. All of the 
quantities are uniform in the y direction. 

(b) The model is a boundary layer model in which dynamic processes 
are neglected. 

(c) Vertical turbulent transfer is parameterized in terms of turbu- 
lent exchange coefficients. 

(d) In the absence of fog, radiative flux divergence is neglected. 

(e) No microphysical calculations are carried out. All influences 
of the fog microphysics in the model are parameterized in terms of the pre- 
dicted liquid water content. 

(f) Supersaturated water vapor condenses instantaneously until 
saturation is achieved. Liquid water in an unsaturated region evaporates 
instantaneously until saturation is achieved or the liquid water is exhausted. 
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2 . 2.2 


Equations 


• List of Symbols 

In order to avoid lengthy explanations in the text, a list of the 
most important symbols employed will be given first: 


T,0 

r 

^s 

w 

u 

u 

g 

V 

z 

k 

X 

i 

t 

n 

K 

m 


P 

C 

F 

R 

a 

P 

L 

V. 


temperature and potential temperature of air 

water vapor mixing ratio 

saturation mixing ratio 

liquid water mixing ratio 

horizontal wind component in x direction 

geostrophic wind 

horizontal wind component in y direction 
height coordinate 

subscript denoting kth vertical grid level 
horizontal coordinate 

subscript denoting ith horizontal grid column 
time 

superscript denoting nth time step 

turbulent exchange coefficient for vertical turbulent 
transfer of momentum 

turbulent exchange coefficient for vertical turbulent 
transfer of heat and moisture 

density of air 

specific heat of air at constant pressure 
net upward flux of infrared radiation 
Stefan- Boltzmann constant 
air pressure 

latent heat of condensation 

mean terminal velocity of fog drops 

spectrally-averaged mass absorption coefficient of fog for 
infrared radiation (cm^ g"^) 

gravitational constant 
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£ 


Coriolis parameter 

u* f r i ct ion vel oc i ty 

H heat flux in constant flux layer 

*3 

L -U C-pT/kgH Monin-Obukhov scaling length 

s p 

roughness height 

k von Karman constant =0.4 


• Major Equations 


The equations employed in the model for the time rate change of 
potential temperature 0, water vapor mixing ratio r, liquid water mixing 
ratio w, wind component u, and wind component v are: 
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3z 
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8t 
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8t 



3z 






CV^w) 


(3) 


3t 


8u 9 3u ^ 

9x ^ 9z ^^m 9z^ 


+ f-v 


(4) 


9v 

9t 


-u — + — (K -r— ) + £. (u -u) 
9x 9z m 9z g 


(5) 


• Saturation Ad j us tment 

The symbol C denotes a source function for condensation or evapora- 
tion. In the actual model, the finite-difference approximations to Eqs . (1) 
through (5) are integrated for a time step, neglecting condensation or evapora- 
tion. Then, the saturation adjustment procedure developed by McDonald (1963) 
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is applied to the new values of 0, v, and w. Taking into account the heating 
of the air by the release of latent heat of condensation, supersaturated water 
vapor at a grid point is converted into liquid water until saturation is 
achieved. Similarly, taking into account the cooling of the air, liquid water 
at a grid point is evaporated into an unsaturated vapor until saturation is 
achieved or the liquid water is exhausted. 


• Turbulent Exchange Coefficients 

As in earlier fog modeling studies (Pilie et al., 1972; Mack et al., 
1972), the most difficult area in this advection fog modeling study was pro- 
viding a realistic description of the vertical turbulent transfer of heat, 
moisture, and momentum over wide ranges of height and stability. In the 
present model, vertical turbulent transfer is parameterized in terms of flux- 
gradient relations using turbulent exchange coefficients. The turbulent ex- 
change coefficients employed in the model are based on empirical flux-gradient 
relations which have been measured in studies of turbulent transfer in the 
lowest 50 m or so of the atmosphere (Businger et al., 1971; Panofsky, 1974), 
where the vertical fluxes of horizontal momentum and heat can be assumed con- 
stant with height, and the Coriolis force can be neglected. The extension of 
the flux-gradient relations obtained for steady-state, horizontal homogeneous 
atmospheric surface layers to the time varying, nonhomogeneous atmosphere 
boundary layers in the advection fog model imposes one of the prinicpal 
limitations to the realism of the model. Nevertheless, in absence of well- 
established higher order closure methods for the prediction of turbulent 
transfer under these conditions, the K- theory approach was adopted for this 
advection fog model to provide a simple, easily understood approximation 
to turbulent transfer of the type employed in earlier fog modeling studies. 


By definition in the constant flux layer (Lumley and Panofsky, 1964), 
the turbulent exchange coefficient for momentum is given by 


2 

u* 


. 9V 
m 3z 


( 6 ) 
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where the friction velocity u* = \/r/p is defined in terms of magnitude t 

9 V 

of the surface stress vector and is the magnitude of the horizontal 

wind shear. Similarly, the turbulent exchange coefficient for heat is 
defined by 


» = -OS'S, f ■ 

90 

where H is the vertical heat flux, and — is the vertical gradient in 

oZ 

potential temperature. Within the constant flux layer, it can be shown that 
the non-dimensional wind shear can be expressed in the form 


isi iX - A 

u* 3z ~ %'-V 


( 8 ) 


*3 

where L = -u C pT/kgH is Monin-Obokhov scaling length. Similarly, the non- 
s p 

dimensional temperature gradient can be expressed in the form 


T* 3z “ 

1 H 

where the scaling temperature is T* = - . 

By definition, then, within the constant flux layer, the exchange 
coefficient for momentum is 


K = 
m 


u*kz 

s 


and exchange coefficient for heat is 


( 10 ) 


^ = 


u*kz 


Cli) 
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In the advection fog model, these turbulent exchange coefficients 
are multiplied by the factor exp[-8fz/u*] , obtained by Shir (1973) in a numeri- 
cal study of a neutral boundary layer, to provide for a reasonable behavior 
of the turbulent exchange coefficients at heights above the constant flux 
layer in the model. 


2 

In the model, non-dimensional wind shear (J> (^ — ) and temperature 

z ^ 

gradient ) relations are employed which are based upon atmospheric 

measurements near the surface under relatively steady-state and horizontally 

homogeneous conditions (Panofsky, 1974). Under stable conditions or 

the relations 


'm 


= 1 + 


=1^ 

s 


( 12 ) 


are used. 


Under unstable conditions or L <0 relations 

s 


"m 



and 



(1 - 15 

s 


(13) 


(14) 


are employed. Under neutral conditions or H = 0, both the stable and unstable 
relations reduce to 


■^h = ‘*’m = ^ (15) 

In applying the advection fog model to study fog formation at sea, 
the roughness height was taken to be 0.01 cm in accordance with wind 
measurements at sea (Hsu, 1974). The first vertical grid level z^ in the 
model was placed at 10 cm height. In order to calculate the friction velocity 
u*, the heat flux H, and the fluxes of momentum, and moisture between surface 
and first grid level, the non-dimensional gradient relations were integrated 
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between and z^ under the assumption that buoyancy has no influence on the 
gradient relations at heights below 10 cm. Because the grid system of the 
model was aligned so that the geostrophic wind u^ is the x direction, the 
horizontal wind at the first vertical grid level z^ has both an x and y com- 
ponent. Assuming that the horizontal wind components u and v are zero at 
the roughness height z^, the following expression is used to compute the 
friction velocity 



where uCz^^) and vfz^^) are the predicted horizontal wind components at the 
first grid level. Similarly, the stress components in the x and y directions 
are 


^ Vu(zp^ -f- v(zp^ . u(z^ 

[ ln(^)]2 

^0 


^ VuCz,)^ + v(z,)^ • V(z,) 

12 

[ln0]^ 

^0 


(17) 


(18) 


Applying a similar reasoning to the potential temperature gradient 
between Zq and z^, the following expression is used to compute the heat flux 
between the surface and the first grid level 



(0(z^) - 0(Zp)) 


(19) 


where 0 (Zq) is taken to be the boundary value on potential temperature. Simi- 
larly, the flux of water vapor between the surface and the first grid level 
is computed from the expression 


10 



( 20 ) 


F = 
r 


k^p 


1 2 
^0 


y; 


(Zj)^ + v(Zj)^ 


Cr(Zj) - t(Zq)) 


where rCz^) is taken to be the boundary value on water vapor mixing ratio at 
the surface. The turbulent flux of liquid water between the surface and the 
first grid level is computed from the expression 


F 

w 


k^p 


^0 


'/u(zp2 + 


v(Zj)' 


W(Zj) 


C21) 


where the boundary value on liquid water mixing ratio at the surface w(Zq) 
is assumed to be zero. 


In the integration of the advection fog model, the quantities u*, 

H, and are evaluated for every grid column at the end of each time step. 

The turbulent exchange coefficients and are, then, computed for every 

grid point in the model from Eqs . (10) and (15), using the friction velocity 

u* and the scaling length determined for the column. In accordance with 

the current practice (Panofsky, 1974), the turbulent exchange coefficients 

for the transfer of water vapor and liquid water are assumed equal to turbulent 

exchange coefficient for heat Because of the implicit integration procedure 

used to integrate the vertical transfer terms in the finite-difference equations 

of the model, the fluxes H, F , F , t , and t are converted into an effective 
^ ' w r^ x^ y 

turbulent exchange coefficient 


K = K, L. 

m C*^) h C-yj z - 

[In(T^)] 

which is used in the finite-difference equations to compute the turbulent 
fluxes between surface and first vertical grid level z^. 


/ 


u(z^)^ + v(z^)^ 


( 22 ) 
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Making use of Eq. 
(Panofsky, 1974) 


(6) to eliminate u*, and using the relation 






m 


( 23 ) 


2 

to eliminate t — , Eqs . (10) through (15) can be converted into expressions for 
^s 

turbulent exchange coefficients K and in terms of the velocity gradient 

and the gradient Richardson number 

d Z d Z 


Ri = 


0 Bz 


.Bu. 2 


( 11)2 


(24) 


Turbulent exchange coefficients of this type were initially utilized in this 
numerical modeling study of advection fog. A dependence of the turbulent 
exchange coefficients on the predicted local gradients of potential tempera- 
ture and wind was believed to permit a more realistic description of the for- 
mation and development of advection fog than a dependence only on the predicted 
surface fluxes, even though both forms were derived from the same constant flux 
layer relationships. Unfortunately, computational instabilities were encoun- 
tered in the integrations of the two-dimensional advection fog model with the 
turbulent exchange coefficients depending on the local wind shear and Richardson 
number, which could not be eliminated within the cost and time constraints of 
the program. Therefore, the dependence of the turbulent exchange coefficients 
on surface fluxes discussed above was adopted in the advection fog model to 
achieve computationally stable integrations. 

• Radiation 


The treatment of radiation in the present model is designed to cap- 
ture the essence of physical processes while avoiding detailed radiative trans- 

BR 

fer calculations. In the absence of fog, the radiative flux divergence — 

oZ 

in Eq. (1) is assumed to be everywhere zero. In the presence of fog, a 
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radiative flux divergence representing the absorption and radiation of infrared 

radiation by fog drops is introduced into Eq. (1). An approximation expression 

3R 

relating the radiative flux divergence to the vertical distribution of the 
liquid water mixing ratio w(z) in the model is obtained by employing a spec- 
trally-averaged mass absorption coefficient to represent the influences of 
fog on radiative transfer. For simplicity, the fog is assumed to be at the 
absolute temperature of the surface, thereby neglecting the influences of 
temperature gradients on the radiative transfer in the model fog. 

It is assumed that a flux radiation equal to the black body radia- 
tion from the surface is transferred upward through the fog or 

Under clear sky conditions, it is also assumed that a downward flux 
of back radiation from the atmosphere equal to approximately 75% of the black 
body radiation from the surface is incident on the fog top at height 
(Sutton, 1953) or 


(26) 

This downward flux of radiation on the fog top from the atmosphere is deficient 
in energy in the water vapor window portion of the spectrum between approxi- 
mately 8.5 and 12 pm wavelength. 


Inside the fog, the downward flux of radiation is determined by 
absorption and radiation of energy by the fog. At height z in the fog, the 
downward flux of radiation is 

z? 

R^(z) = .75aT^'^exp[-1.6^k^(z')pw(z-')dz'] 

+ f 1.6 k^(?) pw(c)exp[-1.6 k^(z')pw(z^)dz"]d? 


Zt 

0 " -^SexpC-l.e / k^(z-')pw(z')dz"]J 


(27) 
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The first term on the right represents the downward flux of radia- 
tion incident on the fog top as attenuated through absorption by the fog in 
the distance from to z. In the second term on the right, the downward 
radiation reaching height z due to radiation from an element of fog of thick- 
ness dc at height C is integrated from height z to the fog top at height z^. 

Here the product pw(z) of the density of air times the liquid water mixing 
ratio is the liquid water content of the model fog at height z. The influence 
of the angular dependence of the radiation field has been approximated by 
increasing the optical path lengths by the diffusivity factor 1.6 (Goody, 1964). 

The net radiative flux at height z in the fog becomes 

T 

k^(z")pwCz")dz"] (28) 

Differentiating with respect to z, we obtain the expression used in the model 
for the radiative flux divergence produced by fog 


RdCz) = 


9z 


(z) 


.25aT^ 1.6k^(z)pw(z)exp[- 


k^(z")pw(z")dz"] 


(29) 


To determine the form for the mass absorption coefficient k^ employed 

in the model, we must consider the absorption cross sections of fog drops. For 

fog drops of radius greater than approximately 4.5 ym, Zdunkowski and Nielsen 

(1969) have shown that the spectrally-averaged absorption cross section of a fog 

2 

drop for black body radiation is within 20% of irr , where r is the radius of 
the drop. For a fog drop concentration of N drops per unit volume, the absorp- 
tion per unit length in the model fog is assumed to be 


1.6 k^pw 


1.6 Nirr 


2 


(30) 


where r is the mean radius of the fog drops. 
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Solving for the mean fog drop radius r from the relation 


Np^ 


47T 3 


pw 


we obtain 


r = 


f 3pw 


03 


(313 


(323 


Substituting for r in Eq. (303, we obtain the expression 


1 i i 

1.6 k^w = i.6(7rN3^(^3^pw3^ 


(333 


for the absorption per unit path length in the fog. 


3 

Assuming a constant drop concentration of N = 50 drops per cm , 
reasonably representative of coastal advection fogs (Mack et al., 1975), 

Eq. (33) reduces to expression 

2 

1.6 k^^pw = 7.1Cpw3\m"^ (343 


for the absorption per cm employed in the fog model as a function of the pre- 

3 

dieted liquid water content pw in grams per cm . To calculate the radiative 
9R 

flux divergence in the potential temperature equation of t^e model, the 
quantity 1.6 k (z) w(z) in Eq. (29) is replaced by 7.1(pw(z))^, and the verti- 
cal integral is evaluated numerically. While this procedure is a crude approxi- 
mation to radiative transfer in fogs (see Korb and Zdunkowski (1970) for 
detailed radiative transfer calculations in fog) , it should provide a roughly 
quantitative simulation of "the alteration of the vertical distribution of net 
radiation by fog and the accompanying radiative cooling of the fog. It has 
the advantage of imposing no significant computational requirements. 
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• Terminal Velocity of Fog Drops 


In the model, the sedimentation of the fog drops is simulated through 
the mean terminal velocity in Eq. (3) , Under the assumption that the total 
drop concentration in the model fogs is constant, is a function of the local 
liquid water mixing ratio which allows to remain small until the predicted 
liquid water contents approach values observed in well-developed fogs. 


The liquid water mixing ratio w can be written 


w 




r 


3 


(35) 


where N is the number of drops per unit volume, r is the mean volume radius 
of the drop-size distribution, p j is the density of liquid water, and p is 
the air density. Eq. (27) can be solved for r and the result substituted in 
the Stokes relationship (Fletcher, 1966) 


= 1.2 X 10^r^(cgs units) 


(36) 


for terminal velocity of water droplets under 20 ym in radius, 
expression is 

2 _ 

= 5.3 X 


The resulting 


(37) 


where N is number of drops cm 


-3 


_3 

Assuming a constant drop concentration N = 50 cm , reasonably repre- 
sentative of coastal advection fogs (Mack et al., 1975), Eq. (37) reduces to 
the expression 


\ = 


2 3 

4 X 10 w cm/sec 


(38) 


-4 

employed in the model. For a liquid water mixing ratio w = 2.44 x 10 corres- 

_3 

ponding to liquid water content of approximately 0.30 g m , Eq. (38) yields 
= 1.5 cm sec"^. 
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Fairly realistic results have been achieved in this study and earlier 
fog modeling studies (Pilie et al., 1972; Mack et al., 1972) using this para- 
meterization of the influences of fog drop sedimentation on fog development. 
However, since the drop sedimentation term has been demonstrated to have a. 
large influence on predicted fog development in the model, careful considera- 
tion should be given to a better description of fog drop sedimentation in 
future fog modeling studies. The description of fog drop sedimentation in the 
model might be significantly improved by carrying out even-highly parameterized 
microphysical calculations in the model. 

• Horizontal Advection 

The horizontal advection terms in the differential Eqs . (1) through 
(5) are approximated in the model by upstream differences which have the form 


M ~ Qj-l,k 

Bx i,k ^ \ \ ^ 


C39) 


for u. , >0 . 
i,k- 

When combined with forward time differencing, this finite-difference 
scheme is computationally stable when 


u 


Ax 


:< 1 


(40) 


Thus, the maximum value of u and the minimum value of the horizontal grid 
length Ax determine the maximum value of the time step At for which the inte- 
gration will be stable. If u = 6.0 m sec and Ax . = 1000 m, then At 

* max min max 

»167 sec. 


The upstream differencing scheme has well-known pseudo-diffusive 
properties (Molenkamp, 1968), which can be expressed in terms of a horizontal 
pseudo-diffusion coefficient 
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C41) 


No unique value can be assigned Vj^ in the model because u varies in the verti- 
cal and Ax varies in the horizontal in the expanding portions of the horizontal 

2 2 -1 

grid. Above a few meters height, however, generally exceeds 10 m sec in 
the numerical experiments to be discussed. This results in the rapid damping 
of short wavelength disturbances in the horizontal. While this property of 
the finite difference scheme would be particularly unattractive in studies of 
transient solutions, it is not without advantage in the study of forced, steady- 
state solutions, particularly if the present model is generalized in the future 
by the incorporation of prognostic fluid dynamic equations. 

2.2.3 Initial Conditions 


At each vertical grid level in the two-dimensional advection fog 
model, the initial values of the prognostic variables 0, r, w, u, and v are 
normally uniform in the horizontal. The model permits initialization with 
any vertical distributions of the prognostic variables. In the fog formation 
cases to be discussed, the model was initialized with vertical distributions 
generated by allowing a one- dimensional (in z) version of the model to evolve 
to approximately steady-state above an ocean surface about 12 hours after 
being initialized with some simple theoretical vertical distributions of the 
prognostic variables. These initial vertical distributions generated by the 
one-dimensional model can be thought of as corresponding to an air mass which 
has been conditioned by a significant trajectory over a uniform ocean surface. 

For testing purposes and simple fog formation simulations, the advec- 
tion fog model is easily initialized with certain idealized vertical distribu- 
tions of the prognostic variables. The initial distribution of potential 
temperature can be either adiabatic (0 = constant) or correspond to an iso- 
thermal temperature distribution. The water vapor mixing ratio can be 
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initialized with a uniform distribution with height. The liquid water mixing 
ratio w can be initialized to be zero everywhere. The u component of the 
horizontal wind can be initialized to have a logarithmic height dependence 
of the form ln (^ ^ component of the horizontal wind can be ini- 

tialized to have a logarithmic height dependence near the surface, but to 
approach zero at the upper boundary. 


2.2.4 


Boundary Conditions 


In the numerical experiments to be discussed, the water vapor mixing 
ratio at the surface was maintained equal to the saturation mixing ratio at 
the temperature of the surface. This boundary condition was applied to study 
fog formation and development over an ocean surface. The liquid water mixing 
ratio at the surface was maintained equal to zero on the assumption the ocean 
surface represents a sink for atmospheric liquid water. At the surface, both 
the u and v components of the horizontal wind were set to zero. 

The basic objective in fog formation and other numerical experiments 
with the advection fog model is to study the response of the model to discon- 
tinuities in the surface temperature. Since the initial condition is a uniform 
surface temperature, the surface temperature discontinuities are introduced 
gradually over the first 30 minutes of the numerical experiments to avoid 
sudden pulsing of the model at the start of an integration. The model allows 
any structure of the surface temperature to be introduced. 

The upper boundary of the grid system at approximately one kilometer 
height is assumed to be an undisturbed level. There, the prognostic variables 
are maintained equal to their initial values throughout a numerical experiment. 
The value of the u component of the horizontal wind at the upper boundary is 
assumed to be the geostrophic wind u^ in the equation of motion for the v 
component. The v component of the horizontal wind at the upper boundary is 
correspondingly maintained equal to zero. 
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At the upwind boundary of the grid system, the normal boundary con- 
dition is to maintain the prognostic variables equal to their initial values 
during an integration. At the downwind boundary, the values of the prognostic 
variables are continually adjusted to be equal to the computed values in the 
adjacent upwind column of the grid system. While there is no mechanism for 
upwind propagation of information in the present model, this procedure was 
adopted with an eye to eliminating reflections of gravity waves at the down- 
wind boundary if the model is generalized in the future by incorporation of 
dynamic process. 

2.2.5 Computational Procedure 

• Grid System 

In the vertical grid employed in the model, the spacing between 
adjacent grid levels expands upward from the surface. The expanding grid 
system provides high resolution near the surface where the variables of the 
model change rapidly with height and removes the upper boundary from the 
region of primary change, without requiring a prohibitively large number of 
grid levels. In the numerical experiments to be discussed, 43 grid levels 
have been employed with the grid spacing expanding upward by a factor of 

1.2 per level from an initial spacing of 10 cm between the first grid level 

in the atmosphere and the surface. The upper boundary was 205 m above the 
next highest grid level and 1060 m above the surface. 

In the horizontal grid system, the grid expands both upwind and 

downwind about an interior region with uniform grid spacing. Within the 
limitation of a maximum of 40 horizontal grid columns, the model permits 
great flexibility in specifying the horizontal grid parameters for different 
applications . 
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• Implicit Integration 


Since the vertical grid has very fine grid spacing near the surface, 
it was necessary to adopt an implicit treatment of the vertical diffusion 
terms in the partial differential equations of the model (Eqs . (1) - (5)) 
in order to obtain computationally stable solutions using reasonably large 
time steps. 


Omitting the symbolic source term C for condensation, the finite- 
difference equations employed in the model can be written in the form 


- 0" 

At 




n 60“ ^ 6 , n 


2 

1 ^ 1000 .^ 
pC P ^ Bz 
P 


(42) 


- r" _ n ^ 6^"' 

At 6x 6z ^ n 6z 


(43) 


n+1 n 
w - w 

At 


^ ^ ^ n+1 6(VV^^) 

n 6w .6 ow ^ . t 


__ _ 

^ ^ ^ ^ ^ ^ 


6z 


(44) 


n+1 n .n . .n+1 

u - u n 6u 6 .^n 6u . jr n 

-u -^ + ^ (K — ) + f-v 

6x 6z m 6z 


At 


(45) 


n+1 

j 

At 


- u^ n 6u^ 6 . n 6v^^^. 

= “U ^ + — (K^ — ) 


6z ^ m 6z 


+ f 


CUg - u“) 


(46) 


where n and n+1 denote values known at the end of successive time steps. 

3R 

Thus, the horizontal advection terms, the radiative flux divergence — , 

d Z 

the exchange coefficients and Kj^, the mean terminal velocity of the fog 
drops V^, and the Coriolis terms in the horizontal wind equations are computed 
from the known values of the prognostic variables at the end of the previous 
time step n. The prognostic variables in the vertical diffusion terms and 
w in the drop sedimentation term are the new values to be determined at 
the end of time step n+1. 
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The method of solution of these implicit finite-difference equations 
of the advection fog model by the technique of Richtmyer (1957) is outlined 
in the discussion of an earlier version of the model by Mack et al. (1972). 

The implicit integration method and the computational sequence employed in 
the advection fog model are fundamentally unchanged by the incorporation of 
the prognostic equations for the horizontal wind and other improvements to 
the model accomplished during this investigation. The reader is, therefore 
referred to the report by Mack et al . (1972) for a discussion of the numerical 

integration method used in the advection fog model. 

2.3 APPLICATION OF THE MODEL 

2.3.1 Introduction 


Numerical modeling studies are at least two-part projects. One is 
the task of developing a reasonable mathematical description of the physics 
of the problem and a stable, accurate numerical representation of that mathe- 
matical description. The other consists of applying the model to the phenomena 
under study and interpreting the results. Of fundamental importance to both 
tasks is a concept of how the phenomenon forms and develops. In addition to 
describing the physics involved, it is necessary to provide a combination of 
realistic initial and boundary conditions from which the phenomenon can evolve. 
These concepts develop through combined field investigations and theoretical 
modeling studies. 

The original development and testing of the Calspan advection fog 
model (Mack et al. , 1972) was carried out under woeful lack of observational 
information. Since 1972 Calspan has conducted observational studies of marine 
fog off the coast of California, which can provide guidance for the develop- 
ment and application of numerical fog models. These studies have demonstrated 
that marine fog forms and develops by a variety of mechanisms (Mack et al . , 
1975). One observed mechanism, stratus lowering, was simulated in a numerical 
modeling study of marine fog by Barker (1973). In the present study, the 
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advection fog model has been applied to simulate the formation of marine fog 
by boundary layer exchange processes, which result when nearly saturated air 
is advected over a change in water temperature. In the marine fog program, 
local fogs have been observed to form when nearly saturated air was advected 
over warmer water. Since the representation of vertical turbulent exchange 
in the model precluded the proper simulation of fog formation by the advection 
over warmer water, the formation of marine fog by advection of nearly saturated 
air over colder water was investigated with the advection fog model. Where 
possible, the marine fog observations were used to guide the investigation 
and delineate important areas for future improvement and application. 

The present study is directed at investigating marine fog formation 
in response to water temperature discontinuities. In order to assure that 
the simulations with the advection fog model would respond soley to this type 
of surface forcing, upwind vertical distributions of the model variables were 
employed which were very nearly in steady-state equilibrium with each other 
and the upwind surface temperature. These upwind distributions were generated 
by initializing a one-dimensional version of the model with some idealized 
vertical distributions of the variables, and by running the model until the 
vertical distributions evolved to a quasi-steady state, compatible with the 
boundary condition aloft and the simulated ocean boundary conditions at the 
surface. In this manner, upwind distributions were employed in the advection 
fog model which were representative of air which had been conditioned by a 
significant trajectory over a uniform ocean surface. These distributions 
were advected over decreases in surface temperature in the advection fog. model 
to investigate the formation of marine fog by that process. 

2.3.2 Fog Formation by Mixing 

To simulate marine fog formation as a result of water temperature 
discontinuities, it is assumed in the model that the boundary condition on 
the water vapor mixing ratio at the surface is saturation at the surface 
temperature. Downwind from a decrease in surface temperature, saturated air 
at the new surface temperature is mixed with warmer, nearly saturated air 
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advecting over the surface temperature decrease. It is this mixing process 
which is responsible for initial fog formation in the advection fog model. 

The water temperature decrease that can produce a fog by the mixing process 
can be determined from the upwind conditions with the aid of a saturation 
mixing ratio vs. temperature curve in a manner first described by Taylor 
(1917) and illustrated in Figure 1. 

In reference to Figure 1, assume that the upwind air advecting over 
a decrease in water temperature is at temperature T and at dew point T^^. Let 
this air mix with saturated air at temperature T^^ over the surface downwind 
of temperature decrease. If heat and water vapor are mixing in an identical 
manner, as in the advection fog model where the turbulent exchange coefficients 
for heat and water vapor are assumed to be equal, all points lie on a straight 
line between the unmodified upwind air and the saturated surface air downwind 
of the temperature decrease. In this case, the mixing process does not produce 
supersaturation and fog cannot form. If the downwind water temperature is T^, 
however, significant supersaturation is produced by the mixing process and fog 
will form. Supersaturation is produced in the mixing zone at all downwind 
water temperatures colder than T^, the temperature at which a straight line 
passing through the point defined by the upwind conditions is tangent to the 
saturation mixing ratio curve. Any water temperature decreases smaller than 
T-T^ will result in no fog formation by the mixing mechanism. It is possible, 
therefore, from such arguments to estimate the minimum change in a surface 
temperature required to produce a fog in the advection fog model for given 
upwind conditions, i.e . , temperature and dew point. 

2.3.3 Upwind Boundary Conditions for Two-Dimensional Simulations 

Consideration of fog formation by mixing indicates that the dew 
point depression upwind of a sea surface temperature discontinuity determines 
whether fog will form as a result of advection over the discontinuity in sur- 
face temperature. For a pure water surface, it is predicted from mixing 
argmnents that the dew point depression, upwind of an 8^C decrease in surface 
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TEMPERATURE 


Figure 1 SCHEMATIC CURVE OF SATURATION MIXING RATIO VS TEMPERATURE 
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temperature, must be O.S^C or less for fog formation. It is apparent that 
careful consideration of upwind dew point depressions is required in order 
to obtain reasonable results in two-dimensional numerical simulations of 
marine fog formation by cooling. 

Ideally, upwind boundary conditions should be vertical profiles of 
temperature and moisture which were observed in air which later experienced 
fog formation. Such observations, in particular those of moisture, do not 
exist. The approach used in this study was to generate upwind boundary con- 
ditions by running the numerical model in a one-dimensional mode to simulate 
the long overwater trajectory. There were two advantages to this method. 

One, the upwind boundary conditions for the two-dimensional model were in 
mutual adjustment throughout the model depth, i.e., the profiles of the prog- 
nostic variables were in balance with the profile of the turbulent exchange 
coefficient. Two, the moisture profile was generated from the surface boundary 
condition of saturation. 

Running the model in a one-dimensional form still required initial 
conditions. In July 1973 Calspan took aircraft soundings 50 n mi out to sea 
northwest of Vandenberg AFB in California. These observations were taken in 
a 400-600 m thick boundary layer with stratus present in the upper portion. 
Below the base of the stratus cloud, the temperature profile was adiabatic 
indicating well-mixed conditions. The moisture profile also indicated well- 
mixed conditions. These observations suggested that a well-mixed layer might 
be a reasonable initial condition for the one-dimensional simulations. 

When one-dimensional simulations were initialized with a uniform 
potential temperature equal to the surface temperature, the temperature dis- 
tribution did not change with time. The water vapor mixing ratio was initial- 
ized to constant values near the surface corresponding to a few degrees Celsius 
dew point depression. Aloft the mixing ratio was decreased to avoid initial 
supersaturation. As a result of the saturation boundary condition on the mix- 
ing ratio at the surface, evaporation occurred and the mixing ratio increased 
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with time in the course o£ one-dimensional simulations. Because of the rapid 
decrease of temperature with height, stratus formation eventually occurred 
unless initial mixing ratios aloft were very small. 

Although stratus formation and the lowering of stratus to the surface 
under the influences of radiation cooling is an important mechanism of marine 
fog formation off the California coast (Mack et al,, 1975), the objective of 
the present investigation was to examine the triggering of marine fog formation 
by advection of nearly saturated air over water temperature discontinuities. 

It was concluded that, while the well-mixed initial condition on temperature 
was conducive to stratus formation in the one- dimensional simulations, it was 
not conducive to the production of nearly saturated conditions near the surface. 

To obtain smaller dew point depressions near the surface and, at the 
same time, avoid stratus formation, one-dimensional simulations were initialized 
with an isothermal temperature profile. However, because of an initial imbalance 
between the isothermal profile and the profile of turbulent exchange coefficient, 
the eddy heat flux showed both convergence and divergence. With convergence 
in the lower levels and divergence around 100 m, the isothermal profile evolved 
to nearly an adiabatic profile from the surface up to 100 m. Above 100 m the 
profile remained isothermal at the original temperature. One of these one- 
dimensional simulations produced a dew point depression of 0.323^C at the sur- 
face with a second minimum of 0.2^C at the top of the adiabatic layer. 

When these steady-state conditions were used to initialize a two- 
dimensional simulation, numerical stability problems appeared with the Richard- 
son number formalism for the turbulent exchange coefficients. After much work 
to eliminate this instability proved fruitless, the Richardson number formalism 
was abandoned and a constant flux approach was adopted which used surface heat 
flux to specify the stability. With the constant-flux formalism, small values 
of the turbulent exchange coefficients developed in the model above fog formed 
by cooling the air from below. Consequently, study of the complete develop- 
ment of this type of fog to large depths was no longer feasible. However, 
both fog formation in the surface layer and the early stages of fog develop- 
ment could be studied. 
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2 . 3.4 


Fog Formation by Cooling from the Surface 


Starting from an isothermal profile, a one -dimensional simulation, 
which used the surface heat flux foarmalism for the turbulent exchange coeffi- 
cients, produced a slightly lapsed temperature profile with approximately 
0.4^C difference between the surface and 100 ,m. The dew point depression 
was 0.5°C near the surface and increased to l.O^C at 300 m. The saturation 
mixing ratio vs. temperature curve indicated that a temperature decrease of 
at least 8^C would be needed to form fog with the O.S^C dew point depression. 

In a two-dimensional simulation, these upwind boundary conditions were advected 
over a surface which was 8^C colder than the upwind surface temperature. The 
horizontal dimension of the colder water was 20 km. At the wind speeds 
involved in the lowest levels (~1.5 m/sec), a four-hour simulation was needed 
to set up a steady state throughout the 20 km horizontal extent of the model. 

At the extreme downwind location, a 9 m deep fog formed with an 

3 

average liquid water content (LWC) of 0.15 g/m . Due to the large decrease 
in surface temperature, the air temperature profile was strongly inverted, 
with 40 m being the height at which no cooling had taken place. Between the 
surface and 40 m, the temperature difference was 6.6^C producing a mean. gradi- 
ent of 0.165^C/m through the layer. This gradient was an order of magnitude 
greater than that in the inversion observed by Taylor (1917) in the low levels of 
fog off Newfoundland. In an attempt to bring the strength of the inversion 
produced in the model fog closer to the observed inversion, a smaller tempera- 
ture decrease of -4°C was used. From the mixing argument for a -4°C temperature 
change, a surface dew point depression of 0.25^C would be needed as the upwind 
condition. 


To obtain very small dew point depressions in the one-dimensional 
model, u* was decreased to 5 cm/sec to reduce the heating resulting from the 
initial imbalance associated with the isothermal profile. Without this reduc- 
tion, the model could not produce dew point depressions under 0.3^C. The 
temperature profile shown in Figure 2 was obtained from the one- dimensional 
model in this mode; the corresponding dew point depression was 0.167^C at the 
surface . 
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Figure 2 STEADY STATE TEMPERATURE PROFILE GENERATED BY 
ONE-DIMENSIONAL SIMULATION 


29 


Because of the smaller wind speeds at the lower levels (0.75 m/sec), 
only 10 km was needed in the horizontal to encompass the steady-state condi- 
tions after four hours.' The steady-state fog at 10 km was only five meters 
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deep with LWC equal to 0.1 g/m confined to below 1 m. The cooling was con- 
fined to 20 m and below. The temperature difference between the surface and 
20 m was 3.5°C producing a gradient of O.iys'^C/m, the same magnitude as in the 
-8^C case. Although the cooling in the boundary layer was proportional to 
the change in the surface temperature, the cooling was confined to a shallower 
layer, so that the temperature gradient in the inversion remained essentially 
the same. The discrepancy between the model results, and Taylor’s measurements 
may lie in the measurements. Modern, more detailed measurements are scheduled 
to be taken in fog situations off Newfoundland in the summer of 1975. These 
measurements should provide better observations for comparison with the model 
results. 


In summary, the model results show heat loss to a cooler ocean sur- 
face at saturation can produce fog formation, but because of the real boundary 
condition of saturation over a salty ocean, an initial dew point depression 
of about 0.25°C requires a 6.5^C change in the surface water temperature. 
Measurements off the California coast have shown that changes of this magni- 
tude in the surface water temperature are not required where fog actually 
forms. Apparently some process, other than loss of heat to the sea surface, 
is responsible for formation of fog under these circumstances. However, it 
is likely that such large temperature discontinuities may exist on the north 
side of the Gulf Stream, an area scheduled for field investigation in the 
summer of 1975. If fog is found to form there under large decreases in sea 
surface temperature, determination of whether the fog forms by cooling will 
require accurate measurement of humidities near saturation. 

2.3.5 Effect of Drop Sedimentation on Fog Growth 

Although the model fogs formed above did not grow to large depths, 
some insight into the effects of drop sedimentation and turbulent intensity 
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on fog growth can be obtained by comparison of the two simulations. A compari- 
son of the two steady-state fogs is shown in Table I. 



Table I. 



COMPARISON 

OF STEADY-STATE FOG PROPERTIES 
SURFACE TEMPERATURE CHANGES 

FOR -4°C and -8°C 


Fog Top Height 
Cm) 

Turbulent Exchange Coefficient 
at Fog Top (cm^/sec) 

Characteristic 
Water Content 

Liquid 

(g/m3) 

-4°C 4.85 

37.7 

<0.1 


1 

00 

o 

n 

00 

267 

>0.1 



For the fog with the higher top, the turbulent exchange coefficient 
is about six times larger than that with the lower fog. Larger K’s were pre- 
sent at the fog top not only at 20 km, but at all locations between 20 km and 
the upwind boundary. The higher K's in the -8°C case arose from higher K*s 
in the upwind condition. Qualitatively, it might be expected that in the 
case with the larger K*s the fog might grow to higher heights than it did. 
However, the sedimentation is larger in the -8°C case. In the -8°C case, the 
fog drops fall at approximately 1 cm/sec, while in the other case they fall 
at approximately 0.5 cm/sec. Thus, in the -8^C case, the promoting and inhi- 
biting processes for fog growth are both larger, with the net effect being 
only a slightly higher fog. 

Further insight into the control of fog development by sedimentation 
was obtained by running a simulation in which the only change was elimination 
of the sedimentation. This type of simulation was run for the conditions of 
the -4°C case, and the fog grew to 70 m in depth at the 10 km location. How- 
ever, this growth was not entirely due to absence of sedimentation. The tur- 
bulent exchange coefficients were also larger. In the model, the strength 
of the K’s is related to the liquid water content through the radiative cooling. 
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If larger liquid water is present over greater depths, less stable temperature 
profiles are produced and, hence, strong turbulence occurs. However, the K*s 
in the simulation without sedimentation were the same order of magnitude as 
those in the -8°C case, in which fog grew only to 9 m. Therefore, the major 
effect on fog development was the absence of sedimentation. Since sedimenta- 
tion has such a profound effect on fog development, the fall velocity of drops 
must be modeled carefully in any numerical simuations designed to study fog 
growth . 

2.3.6 Growth of Shallow Fog Over Warmer Water 


The modeling work has shown that fog formation is possible over the 
ocean by cooling in a shallow layer near the ocean surface. However, the 
large discontinuities in surface water temperature that are required have not 
been found in the observations of fog occurrence over the open ocean west of 
California. Analysis of that field data (Mack et al . , 1975) lead to the fol- 
lowing description of fog formation and development by processes other than 
cooling from the surface; 


^'Local fogs have been observed to form in cool, nearly 
saturated air advecting over warmer water. Subsequent tur- 
bulent exchange and enhanced evaporation from the sea surface 
lead to mixing of warm, moist surface air with cool, moist 
air at higher levels and initial condensation in a shallow 
layer. Radiative cooling of this thin layer lifts the inver- 
sion base from the sea surface and further fog development is 
promoted by radiative cooling and enhanced mixing beneath 
the locally induced, low-level inversion." 


The capabilities of the numerical model were utilized to examine the feasibility 
of fog development by these processes. The simulation involved moving a shallow 
fog (formed by previous cooling) over warmer water to determine the effect on 
fog thickness. 

The experiment design was as follows. A 7 m deep fog, which had 
been generated by a water temperature decrease of 4°C, was used as upwind 
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boundary conditions for a two-dimensional simulation. These conditions were 
advected over a water temperature which was 2°C warmer than that which the 
air previously had been exposed to, A control experiment was run in which 
the same boundary conditions were moved over water at the same temperature 
to which the air had been exposed for the previous four hours. 

In the heating case the fog top grew to 18.7 m compared to a height 
of 8.75 m in the continued cooling case. Exposure to heating raised the fog 
top height to a level more than twice that which occurred with continued cool- 
ing. One reason for the difference was the increased turbulent intensity at 
the fog top as the air moved over warmer water. In the heating case, the 
turbulent exchange coefficients at and above the fog top were an order of 
magnitude larger than those in the continued cooling case. Thus, the increased 
turbulence moved the fog top upward by diffusion of both the liquid water and 
the cooling produced by radiation at the fog top. However, an attempt to 
develop a deeper fog by increasing the temperature difference to 3^C led to 
complete evaporation of the fog. As with most parameters involved with fog 
formation, the small change in water temperature increase had a drastic 
impact on the fog behavior. 

The upwind temperature profile and the temperature profiles through 
each of the two fogs are shown in Figure 3. In the cooling case, the tempera- 
ture shows a small decrease uniformly distributed below 20 m, thus maintaining 
the strongly inverted temperature profile. On the other hand, in the warming 
case, temperature increases up through 7 m and decreases from there up through 

30 m. The initial inversion of 3.5°C between the surface and 20 m is reduced 

to a little more than a degree. In this case, the maximum cooling located 
around 14 m is associated with radiational cooling which extends up to the 

fog top at 18.7 m (Figure 3). This result indicates the importance of radia- 

tional cooling to fog development as hypothesized by Mack et al. (1975). 
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LEGEND 


SOLID LINE • UPWIND PROFILE 
DASHED LINE • HEATING CASE 
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Figure 3 VERTICAL PROFILES OF TEMPERATURE AND RADIATIVE COOLING IN FOG 
EXPOSED TO 1) CONTINUED COOLING AND 2) HEATING 
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2.4 


CONCLUSIONS AND RECOMMENDATIONS 


In this study, a two-dimensional numerical model of advection fog 
was improved by incorporation of prognostic equations for the horizontal wind 
and by utilizing both the predicted wind and potential temperature distribu- 
tions to determine the turbulent exchange coefficients in the model. Numerical 
simulations of advection fog formation and development over an ocean surface 
were compared with Calspan*s observations of marine fog off the California 
coast. When air which had been conditioned by simulation of a significant 
trajectory over a uniform ocean surface was subjected to a decrease in surface 
temperature, the advection fog model predicted that a surface temperature 
decrease on the order of 4°C or larger was required for fog formation. Cal- 
span*s observations of marine fog formation in such a conditioned air mass, 
however, indicate that fog formation occurs with surface temperature changes 
of about l^C. In addition, the subsequent vertical development of fog is 
often significantly greater than that predicted by the advection fog model. 

Since the advection fog model is a boundary layer model in which it 
is assumed that turbulent exchange coefficients for heat and water vapor are 
equal, there are at least two possible explanations for the apparent discre- 
pancy between the requirements for fog formation in the advection fog model 
and in the real atmosphere. One possibility is that mesoscale circulations 
induced by advection over surface temperature changes, which are neglected 
in the advection fog model, may promote fog formation in the real atmosphere. 
Calspan*s observations indicate that mesoscale convergence can be an important 
element in the development and persistence of coastal fog (Mack et al., 1975). 

A second possible explanation is that the transfer of heat and water vapor 
between ocean and atmosphere in response to changes in the ocean surface tem- 
perature may proceed at slightly different rates, promoting fog formation 
under certain circumstances. Indeed, Fukuta and Saxena (1973) have argued 
that the greater rate of molecular diffusion of water vapor than of heat in 
the laminar sublayer at the ocean surface aids fog formation over warmer water 
and inhibits fog formation over colder water. In addition, if the turbulent 
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exchange coefficients for the vertical transfer of heat and water vapor are 
slightly different, a similar type of rate dependent (kinetic) effect on fog 
formation would result. 

While the incorporation of dynamic influences into the advection 
fog model is an important objective, there appears to be sufficiently unre- 
solved questions concerning possible kinetic influences on fog formation and 
the influences of turbulent transfer on vertical development of fog to justify 
additional research with the advection fog model in the boundary layer form. 

It is recommended, therefore, that a numerical modeling investigation should 
be carried out of possible kinetic influences on advection fog formation over 
the ocean, using Calspan observations of marine fog formation to guide the 
study. In order to facilitate this investigation and to develop a better 
modeling framework for simulating the vertical development of advection fog 
under the influences of radiative cooling, it is further recommended that a 
more physically realistic closure for turbulent transfer be incorporated into 
the advection fog model. This could take the form of incorporating a prognos- 
tic equation for the turbulent energy and assumed proportionalities between 
the turbulent fluxes and the turbulent energy (Pepper and Lee, 1974). Alter- 
natively, a prognostic equation for the turbulent exchange coefficients them- 
selves could be incorporated into the model (Nee and Kovasnay, 1969) . 
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Section 3 


RECOMMENDATIONS FOR SIMPLIFIED INDIVIDUAL 
ZERO-GRAVITY CLOUD PHYSICS EXPERIMENTS 

3.1 INTRODUCTION 

The following brief discussion is intended to develop recommendations 
to NASA for selected cloud physics experiments applicable to zero-gravity con- 
ditions. As outlined in the work statement, the basic criteria for the selec- 
tion of the experiments are, first, their pertinence to current fog and fog 
modification research and, secondly, their simplicity, best described as 
"carry-on" or "suitcase" experiments. 

Preparations for a complex cloud physics laboratory (CPL) intended 
as a payload for space shuttle sorties have been described recently by Greco 
et al. (1974) and summarized by Greco and Turner (1975). In cooperation with 
the cloud physics community, the 20 most promising experiments have been selected 
and conceptual designs of the experimental setup, including five different cloud 
chambers, have been created (see Table II). By integrating a multitude of test 
devices into the GPL, considerable savings are possible and duplication of 
numerous common components can be avoided. However, since several of the more 
demanding experiments may require substantial additional development, it would 
be advantageous to choose a few comparable tests which can be performed with 
the least elaborate preparation and be developed into simple, individual carry- 
on kits. This would not only yield early scientific results but also provide 
further experience which could be applied in the final design of the GPL. 

Since a great deal of thought had already been invested in drawing 
up the list of 20 most desirable experiments, no additional experiments will 
be taken into consideration here. Rather, we will select from the list of 20 
those that fulfill, in our opinion, the aforementioned criteria of being 
adaptable as carry-on experiments and also have applications to fog physics 
research. 
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Table II. SUMMARY OF EXPERIMENT CLASSES 
(after Greco et al., 1974) 



Experiment Class Number and Title 

Primary 

Chamber* 

Alternate 

Chamber* 

1. 

Condensation Nucleation 

CFD 

E 

2. 

Ice Nucleation 

SDI 

E 

3. 

Ice Multiplication 

SDI 

E 

4. 

Charge Separation 

SDI 

G 

5. 

Ice-Crystal Growth Habits 

SDI 

E 

6. 

Scavenging 

SDI 

G 

7. 

Riming and Aggregation 

SDI 

G 

8. 

Droplet- Ice Cloud Interactions 

SDI 

E 

9, 

Homogeneous Nucleation 

SDI 

E 

10. 

Collision- Induced Freezing 

SDI 

G 

11. 

Saturation Vapor Pressure 

SDI 

E 

12. 

Adiabatic Cloud Expansion 

E 

- 

13. 

Ice Nuclei Memory 

E 

SDI 

14. 

Terrestrial Expansion Chamber Evaluation 

E 

- 

15. 

Condensation Nuclei Memory 

E 

SDL 

16. 

Nuclei Multiplication 

G 

E 

17. 

Drop Collision Breakup 

G 

SDI 

18. 

Coalescence Efficiencies 

G 

SDI 

19. 

Static Diffusion Chamber Evaluation 

SDL 

- 

20. 

Unventilated Droplet Diffusion Coefficients 

SDL 

E 

*CFD 

SDI 

SDL 

= Continuous -flow diffusion 
= Static diffusion, ice 
= Static diffusion, liquid 

E = Expansion 
G = General 
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Our approach to the experiment selection will proceed in the follow- 
ing sequence: First, a list of experiments relevant to fog will be established; 

subsequently, this list will be narrowed down on the basis of zero-g applica- 
bility. Finally, consideration of equipment problems in view of the simplicity 
requirement will further reduce the list of eligible experiments. 

3.2 PROBLEMS IN FOG MICROPHYSICS 

The following brief review of current problems in fog research is 
intended as a guide for selecting experiments relevant to fog. Although a 
later section deals with applicability to zero-g, we will limit the discussion 
at this point to problems in microphysics suitable for experimentation in 
zero-gravity. 

In the initial stages of fog formation, the most important consid- 
eration involves the nucleation and subsequent condensational growth of the 
fog droplets. For a particle of given size and chemical properties, the 
necessary supersaturation for activation as measured with a thermal diffusion 
chamber has been found to differ considerably from what theory predicts; e.g., 
in the case of NaCl, experiments indicate that a 5 times higher supersatura- 
tion is required (Katz and Kocmond, 1973), whereas nominally hydrophobic 
particles have been found to act as condensation nuclei at much lower super- 
saturations than expected (Ruskin and Kocmond, 1971). At the present time 
it is not clear whether the discrepancies can be attributed to surface con- 
taminants on the nuclei, to unaccounted for conditions in the cloud chamber 
or to errors in the particle sizing equipment. 

For a given CCN spectrum and cooling rate of the air, it should be 
possible to accurately predict the ensuing drop size spectrum of the fog. In 
general, however, calculations of drop size distributions in clouds are in 
poor agreement with actual observations. In particular, there have only been 
a very limited number of case studies (see, e.g., Fitzgerald, 1972) where actual 
measurements of drop sizes and nucleus concentrations at cloud base have been 
compared with computed values. Although agreement was encouraging with regard 
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to average drop sizes, discrepancies in the width of the distributions and 
the total concentrations indicate that important problems still exist. Uncer- 
tainties in the measurements point to a need for well-controlled condensation 
experiments which will provide more reliable data for comparisons with exist- 
ing models. 

In the further evolution of fog (past the early condensational stages), 
drop sizes may be reached where coalescence starts to play a role. Although 
this regime is thought to be most important in convective clouds, we have 
encountered numerous situations where shallow fog (i.e., 200-400 m depth) 
produced considerable drizzle only conceivable as a result of droplet growth 
by coalescence. Drizzle has been shown to play an important role in the down- 
ward transport of moisture in the evolution of fog from lowering stratus (Mack 
et al., 1974). Generally, understanding of the coalescence phenomena has pro- 
gressed considerably in recent years, mainly due to a large number of well- 
designed wind tunnel investigations at UCLA. The most difficult area of experi- 
mentation, i.e., that involving relatively small drops, is the one most appli- 
cable to fog; here, further research is needed to compare observations of drop 
coalescence with theory. 

Fog drops not only interact with each other but also with other 
particles. Since, in most cases, only about 10% or less of all particles 
act as CCN, those remaining may be scavenged by droplets. The relevance of 
scavenging to the fog problem can be described as follows: The immediate 

influence of scavenged particles (or gases) on the host droplet is generally 
small, except in the case of surfactants which could retard evaporation of 
the droplet and thus stabilize the fog (see Kocmond et al., 1972). Uptake 
of salt particles or of gases which react in the liquid phase to produce 
salts will cause a slight reduction in water vapor pressure over the drop 
surface in proportion to the increase in concentration of the salt solution. 

The influence, therefore, will be greater as the droplets evaporate. The 
most profound effect of particles and gases being scavenged is in the result- 
ing alteration of the CCN spectrum after complete evaporation of the fog 
droplets has occurred (occasionally referred to as memory effect). Thus, a 
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subsequent fog formed on the CCN spectrum may possess significantly different 
microphysical features. This alteration means enlargement of the previously 
active particles and, in general, a shift of the CCN spectrum to lower critical 
supersaturations ; on the other hand, scavenged surface active substances may 
have a contrary effect by retarding the nucleation of fog droplets. 

A different type of alteration of the CCN spectrum may take place 
when droplets evaporate completely and the residue breaks up into several 
smaller particles. The conditions, other than low humidity, under which such 
a "nuclei multiplication" can occur is not clear at the present time (see 
Podzimek and Saad, 1974). 

Due to the occurence of fogs at subfreezing temperatures, some of 
the additional problems associated with ice nucleation have to be included 
in our considerations of experiments suited to zero-gravity conditions. Since 
most supercooled fogs are stable (i.e., due to a lack of ice nuclei, they 
usually do not precipitate), ice nucleation enters the picture mainly as a 
process employed in fog dispersal. The main difficulty in this type of fog 
abatement scheme is in the economical distribution of suitable ice nuclei 
(i.e., nuclei that are active at relatively high temperatures). By the pre- 
ferred use of homogeneously formed ice crystals (dry ice, liquid propane 
seeding) , it is possible to avoid problems of contact nucleation and ice 
multiplication which are currently among the most debated issues regarding 
the glaciation of convective clouds. 

3.3 RELEVANCE OF EXPERIMENTS TO FOG MICROPHYSICS AND TO ZERO-G 

APPLICATION 

If one examines the list of experiments (Table 1, taken from Greco 
et al., 1974) in light of the above review of problems in fog microphysics, 
the following experiment topics seem to be relevant to this research area 
(in approximate order of importance) : 
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Table III. EXPERIMENTS RELEVANT TO FOG RESEARCH 


1. Condensation Nucleation 

2 . Scavenging 

3. Coalescence Efficiency 

4 . CCN Memory 

5. CCN Multiplication 

6. Static Diffusion Chamber Evaluation 

7. Ice Crystal Growth Habits 

8. Ice Nucleation 

9. Drop let- Ice Cloud Interactions 


From the previous discussion, it is evident that the last three 
experiment classes are not as important to fog research as the first six 
topics. The nine experiments shown in Table III are briefly considered below 
in terms of applicability to zero-gravity tests. 

1. Condensation Nucleation : Since droplet growth at low supersat- 

urations (which are typical of fogs) proceeds very slowly, zero-g conditions 
are especially desirable. As opposed to the situation in the terrestial lab- 
oratory, droplets forming on condensation nuclei can be observed for extended 
times since they will not fall out of the sensitive volume; there are no 
apparent disadvantages introduced by the weightless conditions. 

Scavenging : The advantage of studying scavenging under zero-g 

conditions is two-fold: (1) Elimination of the orthokinetic component will 

help to identify the effects of other scavenging mechanisms and (2) the lack 
of droplet sedimentation will allow for interaction with the particles over 
a time period comparable to actual cloud conditions. 
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3. Coalescence Efficiency ; Droplet coalescence chiefly depends on 
the motion of the drops relative to the air and to each other. In principle, 
therefore, strict zero-g conditions present a disadvantage due to the lack 

of motion; however, by artifically generating a very low gravity situation, 
relatively large drops can be used to simulate the behavior of the more diffi- 
cult to work with small cloud droplets (as proposed by Telford, 1974) . 

4. CCN Memory : Since this effect is essentially a result of scaveng- 

ing, one has to be aware that at zero-g the potentially considerable orthokinetic 
scavenging component will be missing. However, the gravity-free environment 
will permit observation of the same CCN acting before and after uptake of 
material from the surrounding air. 

5. CCN Multiplication : Performance of this experiment would be 

substantially improved if carried out at zero-g, since one could possibly 
observe individual CCN before and after shattering instead of merely record- 
ing a statistical result as is done in terrestrial tests. 

Static Diffusion Chamber Evaluation : In the terrestrial labora- 

tory, one of the drawbacks to making CCN measurements at low supersaturations 
is that the droplets settle out of the sensitive volume before they reach 
detectable sizes. The extent to which this happens at slightly higher super- 
saturations must also be determined. The zero-gravity environment offers a 
unique opportunity to study the fall-out problem and to ‘calibrate’ the static 
diffusion chamber for terrestrial use. 

7. Ice Crystal Growth Habits : Terrestrial limitations are consid- 

erable in this area of research due to the comparatively high fall velocities 
of ice crystals. Zero-g permits prolonged observation without the need to use 
supports for growing ice crystals; but, on the other hand, the lack of venti- 
lation presents a serioias disadvantage (with the possible exception that it 
provides one data point for a growth vs ventilation relationship) . 
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8. Ice Nucleation : Most experiments in this category can be carried 

out successfully in terrestrial chambers due to the generally fast activation 
of ice nuclei. Zero-g would offer an advantage for testing suspected slow 
activating ice nucleating aerosols; however, this application would be some- 
what less important. 

9, Droplet-Ice Cloud Interactions : The interaction between droplets 

and ice crystals consists of vapor transport from droplets to crystals and 
droplet- crystal collisions. The first mechanism is basically covered under 
"ice crystal growth habits", whereas droplet- crystal collisions are largely 

a result of differences in fall velocities. Thus, the gravity-free environ- 
ment would only assist in obtaining data from collisions due to such non-gravity 
effects as, for instance, electrostatic charges. 

The experiments considered so far may be further narrowed by elimina- 
ting the ice-related topics from Table III. Clearly these are important problems; 
however, with respect to warm fog microphysics considerations, they are of less 
relevance . 


3.4 LIMITATIONS OF EQUIPMENT TO SIMPLE CARRY-ON APPLICATION 

The experiments which were selected so far must also be examined 
with respect to necessary and available equipment. In the following paragraphs, 
we therefore review pertinent devices in view of their acceptability for a 
"simple" setup. The main elements of the equipment in question are: 

(1) Test chambers for containment of a specific working environment, (2) 
droplet/particle detectors, and (3) droplet/particle generators. 

1. Test Chambers 

Although some experiments with hydrometeors such as those conducted 
aboard Sky lab IV (concerning the collison of large water drops; Vaughan and 
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Hill, 1974) can be performed without the use of confining walls, it is generally 
necessary to physically define a volume and to control environmental parameters 
in order to maintain suitable conditions. At the present time, we are in prin- 
ciple not restricted in the choice of cloud chambers as indicated by the find- 
ings reported by Greco et al. (1974); however, the basic possibilities are well 
covered by the five chamber types listed in Table II. The chamber with the most 
advantages in principle is the expansion type, since in operation uniform 
temperature changes over the entire volume can be achieved allowing extensive 
relative hiomidity variations from fairly dry to high supers aturat ions . While 
the expansion chamber provides an excellent simulation of atmospheric situations, 
its major drawback from the standpoint of our present task (i.e., simple, carry- 
on equipment) is its complexity in connection with the need to precisely match 
the wall and air temperature. At this stage of technology, this otherwise 
ideal chamber is too complicated to be part of a carry-on experiment. 

Other chambers for consideration include the thermal diffusion chamber, 
whose operating principal (diffusion of water and heat from a warmer upper sur- 
face to a cooler lower one) allows for a humidity regime adjustable between 
saturation and a few percent supersaturation. Although these conditions do 
not provide a true simulation of condensational growth in the atmosphere, they 
offer the possibility to assess CCN activity. The differences between the 
three types of diffusion chambers (listed in Table II) currently employed per- 
tain mainly to temperature range and operating mode. In the continuous flow 
diffusion chamber, CCN are activated while the air sample travels through the 
chamber; while this may have some advantages where recording of the results 
is desired, there are two drawbacks: first, the advantage zero-g offers by 

providing prolonged time for^observing droplets is largely cancelled by the 
motion through the chamber and, secondly, the delicately balanced flow system 
requires careful adjustment of auxiliary flows making the device rather unsuit- 
able for a sin 5 )lified carry-on system. In contrast, both static diffusion 
chambers (liquid surface and ice surface) are operated on the principle that 
measurements are made in situ, i.e., without removing the sample from the area 
of supers aturation, thus affording the benefit of zero-g non- fallout. Whereas 
the chamber with solid (ice) plates may be used for most work on supercooled 
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clouds, the liquid plate chamber is ideal for typical, warm fog-related inves- 
tigations . Both chambers stand out by their simplicity, especially the liquid 
plate version which requires only one temperature control (as opposed to two 
for the ice plate chamber) . 

In contrast to expansion and thermal diffusion chambers in which 
the operating mode provides intrinsically for the formation of hydrometeors, 
the "general" chamber as described by Greco et al . (1974) constitutes a volume 
with temperature-controlled walls; it, therefore, needs additional devices 
for humidification or droplet generation. The apparatus proposed by Telford 
(1974) for the study of droplet-droplet interactions could be considered a 
specialized case of the "general chamber" (basically a fall tube in which 
individual pairs of large droplets (~l/2 mm) are moved under very low gravity 
(-1/1000 g) simulated by rotation of the apparatus). 

From the foregoing discussion, it would appear therefore that the 
static diffusion chamber with liquid plates is best suited for the simplified 
carry-on experiments. 

2. Droplet and Particle Detectors 

The fruitful evaluation of any experiment involving particles and 
droplets requires some knowledge about their concentration and/or size dis- 
tribution. When considering fog microphysics, we are mainly dealing with 
nuclei of 0.01 pm<d<l ym and droplets in the range 0.5 ym<D<30 ym. It is 
current laboratory practice to cover the above particle size range with an 
Electrical Aerosol Analyzer (EAA) and an Optical Particle Counter (OPC; e.g., 
Royco) for nearly real time results or to examine precipitated particle 
samples with scanning or transmission electron microscopes. The latter 
method is very tedious and usually acceptable only for occasional calibra- 
tions ; also, precipitation of particle samples often introduces additional 
problems in the analysis. 
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Both the EAA and the Optical Particle Counter are planned components 
of the complete CPL, but for our present purposes, the following disadvantages 
have to be taken into account: Besides the bulk and weight, the EAA depletes 

aerosol storage due to a relatively high sampling rate (up to 50 liters min*"^) 
and requires higher particle concentrations than generally used in fog- related 
problems, thus creating a need for a dilution system. The OPC only handles 
sizes larger than 0.5 ym which makes it unsuitable as the only particle analyzer. 
An alternate approach to the particle sizing problem will be discussed in the 
paragraph, "particle generation". 

The various methods of assessing cloud droplets also have to be 
tested for their applicability to a carry-on type of experiment. The optical 
methods currently used with diffusion chambers (but applicable to others) are 
direct photography in situ which only records the concentration of drops larger 
than the sensitivity threshold (e.g., 1 ym) , and OPCs which require that the 
drops be moved from their original location to the analyzer; this latter pro- 
cedure is quite critical because the slightest change in relative humidity 
which might occur during transfer will influence the size of the droplets. 

The above-mentioned photographic method has the advantage over OPCs that 
sequential exposures can be taken which may contain information about growth 
rates or which may show the same particle/droplet during several condensation- 
evaporation cycles. Methods of optical probing in situ (i.e., without dis- 
turbing the droplets) by sensing scattered light from a representative sample 
(as opposed to the drop-by-drop procedure of the OPCs) would be most attractive 
in principle, but still suffer from interpretative difficulties. A further possi- 
bility to assess droplet population is to impact the drops on a sensitive film; 
this method, however, also removes the droplets from the test volume and dis- 
turbs conditions inside the chamber. 

If experiments can be designed such that the droplet size does not 
have to be known accurately (as long as it is above a given value) , the photo- 
graphic method is certainly by far the simplest and most straightforward. 
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3. Droplet and particle generation 


As pointed out previously, droplets are easily formed in expansion 
and diffusion chambers provided that suitable nucleating particles are present. 
Thus, the problem is reduced -to producing the desired nuclei. Providing small 
droplets for a '^general" chamber will require a droplet generator of a type 
which yields the particular size spectrum needed for a given application. 

Although the study of natural aerosols as they relate to fogs should 
not be neglected, their complexity makes it difficult to interpret results. 

Also, storage time during transport into orbit would reduce the value of a 
natural aerosol considerably. Among the many techniques for producing arti- 
ficial aerosols, we have to select the ones which utilize simple devices, per- 
form reproducibly and cover the desired size range. In addition, since particle 
analyzers are impractical for our purpose, we will have to work with a particle 
generator which produces an accurately predictable size spectrum and concentra- 
tion--this is practically equivalent to demanding high monodispersity . It is 
desirable to work with soluble and insoluble particles. Several convenient 
methods exist for the production of soluble particles from aqueous solutions 
which are dispersed into droplets and subsequently evaporated. Not all of 
these methods, however, would be suitable for our purpose; e.g., the dipping 
reed device has the disadvantage of not producing a high enough concentration. 
The technique which looks most promising at the present time is the vibrating 
orifice apparatus where a liquid jet of precisely metered flow is broken into 
equal segments by a high frequency vibration of the nozzle (see, e.g.. Wedding 
and Stukel, 1974). The commercially available model can produce several hun- 
dred particles/cc down to less than 0.1 ym diameter. 

Water insoluble particles can be produced with the same device by 
utilizing a non-aqueous solution of the material in question; however, it may 
be difficult to remove the solvent from the atmosphere after the solution 
drops have evaporated. Other methods of generating water insoluble monodis- 
perse aerosols are either too cumbersome for the present purpose (e.g., 

LaMer generator) or not reproducible enough as in the case of simple evapora- 
tion and recondensation. 


50 



3.5 


RECOMMENDATIONS 


From the previous discussion, the following recommendations can be 
made relative to equipment and experiments for use in a shuttle-type spacecraft. 

1, It is recommended that the thermal diffusion chamber be used in 
initial experiments. The chamber should be equipped with a close-up camera 
and illumination should be by laser. This apparatus, in conjunction with a 
vibrating orifice aerosol generator, will allow the following experiments to 
be carried out: 

• Evaluation of static diffusion chamber (SDC) . In view of the 
widespread use of the SDC for measurements of CCN activation spectra through- 
out the world, it would be especially important to perform an experiment to 
validate the SDC. This essential task can be accomplished with relatively 
little effort in the course of initial check out of the equipment at one-g 
and zero-g. It is recommended that tests be performed to determine the extent 
of fall out at low supersaturations and that chamber performance be evaluated 
with several sizes and types of aerosol. 

• Condensation nucleation studies : The size-supersaturation 

relationship for NaCl CCN should be examined carefully in a zero-g environ- 
ment and, more importantly, extended to larger particles than was possible 

at one-g conditions (see Katz and Kocmond, 1973), The same experiment should 
be conducted with CCN of different substances, such as ammonium sulfate, 
which is currently regarded as one of the most abundant CCN materials. The 
obvious benefits of a fallout-free environment and the opportunity for acquir- 
ing significant scientific knowledge make these experiments of highest prior- 
ity. Data generated from these tests will have immediate applications for 
cloud physicists and numerical modelers. 

Testing of a water insoluble aerosol as CCN would be equally impor- 
tant since most natural nuclei are partly composed of some insoluble matter. 
Previously indicated difficulties with solvent removal in connection with 
particle generation could possibly be overcome by use of activated charcoal. 
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• Scavenging, Although a larger chamber would be preferable for 

scavenging studies, it is expected that diffusive losses in the thermal diffu- 
sion chamber would also be acceptable (according to Fuchs, 1964) if particle 
sizes are properly chosen. A tentative experiment could be carried out as 
follows: An aerosol with two distinctly different particle sizes can be 

introduced into the chamber and the supersaturation adjusted for activation 
of the larger particles only. After a certain time, the supersaturation can 
be increased so as to activate the smaller particles, too. By varying obser- 
vation times and growth rates of the droplets grown on large CCN scavenging 
parameters can be deduced. 

2. By slightly modifying the thermal diffusion chamber, it can be 
operated as a "haze chamber"; If a salt solution instead of pure water is 
used on the chamber plates, conditions of slight subsaturation can also be 
achieved. This would make it possible to carry out the following experiment: 

• CCN memory effect . In a procedure similar to the one of the 
proposed scavenging test, droplets formed on the more active particles of a 
two particle type aerosol can be allowed to interact with the non-activated 
particles. Subsequently, the droplets can be evaporated during a brief period 
of sub saturation. Reactivation of the more active particles will show the 
memory effect, 

3. It is recommended that the previously mentioned experiment on 
coalescence proposed by Telford (1974) be considered for development into a 
carry-on package. 

4. Experiments regarding CCN multiplication require that the 
humidity be decreased to very low subsaturation. This is not feasible within 
the diffusion chamber, and the removal of the sample from the chamber for 
drying makes it impossible to observe the same particle/droplet throughout 
the experiment. It is, therefore, recommended that this experiment be post- 
poned until it can be carried out in an expansion chamber. 
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It is obviously impossible, in this initial effort, to develop a 
detailed experimental plan for the carry-on experiments listed above. It is 
recommended, however, that additional consideration be given to fully develop- 
ing the ideas we have presented herein. In particular, the initial testing 
of the static thermal diffusion chamber will be an important first step toward 
the successful completion of the more complicated nucleation studies. The 
further development and testing of a haze chamber may also find special appli- 
cations onboard the CPL. We feel that the experiments suggested here can be 
designed in such a manner as to provide pertinent immediate information for 
later application on the fully developed cloud physics laboratory. 
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APPENDIX A 


COMPUTER PROGRAM DOCUMENTATION FOR ADVECTION FOG MODEL 
A.l INTRODUCTION 

The computer program used in the present study was basically the 
same one used previously and thoroughly documented in Mack et al . (1972). 

Since only minor changes were made in the computer program, the thorough 
documentation is not repeated. Documentation is presented here only for the 
few changes that were made in the computer program. However, for those readers 
who may want to run the program and who may not have access to the previous 
documentation, a description of the input variables and control indices is 
provided which is sufficient for the running of the program. 

Changes from the previous program concern initialization of the two- 
dimensional model from upwind conditions generated by a one-dimensional model 
and computation of the turbulent exchange coefficients. The initialization 
utilizes the read-list capability of the model to which two new control indices 
have been added. The equations for computation of the turbulent exchange coef- 
ficient, which are described in the main text [Eqs . (6)- (15)] are fixed in the 
program and are not modified by input variables or control indices. 

The program operates in cgs units with temperatures in degrees 
Kelvin. All input variables are in cgs units with the exception of the uniform 
spacing between columns (DELX) which is in meters. When liquid water content 
is read in from a list, the units used are grams of water per grams of air. 

The program required 192K bytes of core storage for the execution 
step. The program was run on an IBM 370/168 operating under Release 1.6 VS. 

For one hour of meteorological time, eight horizontal grid columns with 
43 vertical grid levels each, and a time step of two minutes, the execution 
took five seconds of machine time. 
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INITIALIZATION 


A. 2 

Initialization of the wind, temperature, and moisture (vapor and 
liquid) occurs in a loop which is executed two times. The first time through 
all of the columns except the upwind column are initialized to contain iden- 
tical values of the variables. The second time, variables are initialized 
for the upwind boundary (column one) . The initialization procedure is con- 
trolled by input indices where values other than one provide for uniform 
or idealized profiles in the vertical. Index values of one provide for non- 
uniform profiles via the read- list input procedure. 

When the variables are initialized by reading from a list, it is 
possible to initialize potential temperature, water vapor, and liquid water 
from any vertical grid level to the top of the model. This feature allows 
inversions to be specified in the initial profiles with an adiabatic tempera- 
ture and uniform water vapor mixing ratio present below the inversion height. 
The grid level at which the non-uniform values begin are specified by for 
the potential temperature and by ^ for both the water vapor and liquid water. 
In the two-dimensional mode of the model, when a list is provided for initiali- 
zation, KT and KR are set equal to one. 

In the read-list mode of initialization, the type of input is, of 
course, controlled by the format statement. In the program listing presented 
in Section A. 4, a Z-format (hexadecimal) is used since the program is reading 
variables generated by the one-dimensional program. 



SAMPLE SET OF INPUT CARDS 
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A. TIME VARIABLE CONSTANTS, TEMPERATURE AND WIND INPUT ( card ) 
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4444444444444444444444444444 

I 2 J « S I 7 I I II It 12 U H 15 N » II If 70 7171 2J 74 n X 27 n 
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I t I 

77 7 
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I to II 12 11 II IS It 

99939999 

I to It 1 : >1 14 IS II 
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GENERAL PURPOSE - 20 TltlP 
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iiiii 
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]] 34 13 H 
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U34 3SM 
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SJ S4 35 X 
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lllfilllll 
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11111111 
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55S5|5 5 5 5 5 S|S|5S|5 5 5 |5|S5 


fifi fififilfi fifi fi 6 Sfilfi 6 
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2222222 

nMxxnxxiol 
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n nx 

111 


4 444| 
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|«SHI1 II 
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n 71 nl 


CItllt 
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|ti 10 71 17 n 
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11i 
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777 


'4 44 4 
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|n II II u 

7 7 7 7 


1 11 II 1110 

9 9 9 9 


DT time step (seconds) 

OT time interval at which output occurs (seconds) 

ET length of simulation (seconds) 

TIM time interval from t = 0 at which the colder or warmer surface tempera- 
ture reaches its desired value (seconds) 

RF fraction of surface black body radiation used as net upward infrared 
flux through upper boundary 

KW a constant multiplied times liquid water content raised to the two- 
thirds power in equation (34) to obtain infrared absorption per cm 
TP temperature in degrees Kelvin used to initialize the temperature 
profile to isothermal 
ZO roughness length (cm) 

DTEMl change from surface potential temperature in Column 1 needed to obtain 
a uniform surface temperature among the other columns when a uniform 
discontinuity in surface temperature is desired 
UF friction velocity which enters into the quasi-adiabatic velocity pro- 
files when they are used to initialize the velocity field (cm/sec) 
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B. CONTROL INDEX CARD (1 card) 



= -1 PT(I,K) = , isothermal 

IPT = 0 PT(I,K) = PTI, adiabatic 

= 1 PT(I,K) = LIST 


= 0 R(I,K) = RI 

= 1 R(I,K) = LIST 


= 0 WCI,K) = 0 

= 1 W(I,K) = LIST 


ISED = 1 

i 1 


sedimentation 
no sedimentation 


IRAD = 1 

f 1 


radiational cooling 
no radiational cooling 


IF IRSFC = 0 R(I,1) 

= 1 R(I,1) 


R(I,2) 

saturation mixing ratio at surface temperature and 
1000 mb. 
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=0 Quasi-adiabatic u and v profiles 
=1 u and V = LIST 

IDTEM = ° 

= 1 DTEM(I) = DTEMI between ITEML and ITEMR, and equals zero elsewhere 
= -l DTEM(I) = LIST 

-2 -1 

HR = net radiative flux, cal cm min (positive upward) 

COOL = rate of radiative cooling, hr"*^ 

IP = 0 No HR and COOL output provided 
IP ^ 0 HR and COOL output provided 

II = 2 needed in velocity integration when turbulent exchange coefficient 
depends on local vertical gradient of velocity. 

II = 1 when turbulent exchange coefficient does not depend on local conditions. 

KT = vertical grid level to which first potential temperature in list is to be 
assigned. 

KT = 1 when two-dimensional simulation is initialized from list generated by 
one-dimensional simulation. 

KR performs same function for mixing ratio and liquid water content that KT 
does for the potential temperature. 

IZ specifies number of the column for which wind, potential temperature, 
mixing ratio and liquid water will be punched out. 

IZ = 2 in one-dimensional simulations. 


C. GRID SPECIFICATION CARD (1 card) 
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1111 
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ZAL = height of lowest vertical grid point in the atmosphere (cm) 

ZAK = expansion factor in the vertical 

XAI = expansion factor in the horizontal 

IE = number of horizontal columns 

KE = number of vertical levels 

IL = leftmost grid point of unexpanded grid 

IR = rightmost grid point of unexpanded grid 

ITEML = left and right I-limits on non- zero, 

ITEMR uniform DTEM(I) 

DELX = uniform X-spacing (meters) 










D. OUTPUT FORMATS (2 cards) 
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The printed output for the various quantities is in convenient standard units 
as shown in the format heading statements 3000-3800 in the program listing 
in Section A-4. In addition to previous outputs, the outputs of u and v com- 
ponent of the wind, the friction velocity and the eddy moisture flux have 
been added. 







FORTRAN LISTING OF COMPUTER PROGRAM WITH COMMENT CARDS 



c^ vji ^ u> ro 


C 

c 

C TWO-DIMENSIONAL ADVECTION FOG MODEL 

C STEADY STATE MODEL 

C 

C CALSPAN CORP, MARCH 1975 W.J. EADIE AND C.W.ROGERS 
C 

c 

REALMS PT.EPT,FPTfDPT*ATtCTtBT 
REAL*4 KA«INTfKW»LfLM«LHAX»KM 

COMMON PT(40,60)rEPT(60) ,FPT<60) » DPT,AT ,CT, BTf 
e R(40f6G)tW(40t6C)»UC40f60)fT(40t6G)tKA(40f60) t 

7KM(40v60) f 

1NT(40«6C)«CPT(40»60) tHC(4Cf 60) rX(40)«0X(40)«0TEM(40}» 
P(60)»ZA(60) tDZA(60) t ER( 60 «FR (60 ) « EW (60 ) t 

FH(60),EU(60)«FU(60) tZ(60) »DZ( 60 ) t PR (60 ) t L t DEN »CP tG t RA«RW« 
SlGMA,TlMEvDT 9 TlM,UFvRF,KWvZ 0 tCVtUl«DTEMl,ZAK,XAltDELX«CCt 
CHyCItCKtCLf CRfCStUU«FV(60JtEV(60l«V(40960i t F, 

KEyKNylEylNtlLtlRtlSEOflRAOflRSFCtlDTEMflTEMLtnEMRt 1P,II 
DIMENSION FMTHllD.FMTTnU 
C 

C PHYSICAL CONSTANTS 

C 

L«592. 

DEN=1.23E-3 

CP=.240 

GC980.6 

RA=.0686 

RWs.1102 

SIGMA=1.355E-12 

CV=400. 

F:=1.0E-4 

P(l)»1000. 

C 

C DATA INPUT 

C 

10 READC5»1000tEND=400) OT ,0T , ET *TI M ,RFt KW tTP » ZOt DTEMI ,UF 

READ (5 » 1100) IPTf IRRflWtlSEDtlRADvIRSFCflUylOTEMtlPtlltKTvKRi IZ 
READ (5 « 1200) Z AL vZAK »X A I, IE » KE * 1 L «I Rt ITEML t ITEMR * DELX 
REA0(5,1205) FMTH,FMTT 
C 

C DATA LISTING 

C 

WRITE(6»2000) DT »0T , ET , T IM, RF ,KW , TP ,Z0, DTEMI , UF 

WRITE(6«2100) IPTflRRvlWtlSED.IRADtlRSFCtlUtlDTEMflP,!! vKTvKRf 1 
IZ 

MRITE(6«220C) ZALr ZAKt XAI » IE fKE t I Lf IR 9 1 TEHL t ITEMR ,DELX 
C 

C WORKING CONSTANT DEFINITION 

C 

CC=.622*(L**2)/(CP*RA) 

CR=KW*RF*SIGMA/(CP»DEN**.333) 

CI=KW*(0EN*^,667)/2. 

CH=(L**2)/RW 

CK=*.4*UF 

CS=G/4.186E+7 

IN=IE-1 

KNaKE-1 

C 

C VERTICAL GRID SPECIFICATION 
C 
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IA(1|«0 .0 
ZA(2)«ZAL 
DZA(2)«ZAL 
DO 20 K»3fKE 

DZA (K )> ( 1 .4^ZAK ) ♦DZA ( K-1 ) 

20 ZA(K)»ZA(K-1I+DZA{K) 

DO 25 K=2»KE 
Z(K)*(ZACK)+ZA(K-1) )/2« 

25 CONTINUE 

HOR I ZONTAL GR ID SPEC IF 1 CATION 

DO 30 I-ILfIR 

30 Xms(l-IL)4^DELX 
I1«IR^1 

DO 31 I*Il,IE 

31 XU )*XlI-l)^XAI»lXlI-l)-XU-2)) 
IFCIL .EQ. n GO TO 33 
Il==IL-l 

DO 32 I2=1»I1 
I=Il-I2*l 

32 X(I)=X(I+l)-XAI*(X(I+2)“X<I + in 

33 DO 35 1=2, IE 

35 0X(I>=X1I)-X(I-1) 


VARIABLE INITIALIZATION 


READ(5,1300) RI ,PTI 
WRITE(6,2300) RI,PTI 
TP=TP+.l6 
PTI*PTI^,16 
W(l,I)=0.0 

A-UNIFROH WITH HEIGHT 

00 40 1=1, IE 
00 40 K=1,KE 
W(I,K>=0.0 
R(I,K)=RI 
PT(I,K)=PTI 

1F(I ^EQ. n PIKI*1000.*EXP(-G*ZA(KI/(4.186E+7*RA*IP) J 
B-ISOTHERMAL 

IFdPT .EQ. -II PTUfK)=TP*( (1000. /P(K) )♦♦. 286) 

D-QUAS I -ADIABATIC VELOCITY PROFILES 

IFUU .EO. 0) U(I,K) = 2.5*UF*ALOG(UA(K) + ZO)/ZO)*COS(.1745) 

IF( lU .EO. 01 V(I,K)=2.5*UF*AL0G((ZA(K)+Z0|/Z0)*SIN( .1745)*(1-ZA( 
IKI/ZAUE)) 

40 CONTINUE 

INITIALIZE 2-D SIMULATION FROM 1-0 SIMULATION 

IU»2 INITIALIZES UPWIND BOUNDARY III = 1 INITIALIZES REST OF DOMAIN 
DO 500 111=1,2 

IFCIU .EQ. 1) READ(5,4000) ( U( 1 , K ) ,K=X, KE ) 
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IFIUJ .EQ. 1) READ(5t4000) < V( 1 tK ) «K=1,KE) 

KT AND KR PERMIT NQN-UNIFORH INITIALIZATION OF l-DIHENSIONAL 
SIMULATIONS FROM KT AND KR LEVELS TO TOP OF THE MODEL 

IFIIPT .EQ. 1) REAO(S,500D) ( PTC 1 »K) tKsKTrKE) 

IFCIRR .EQ. 1) READ(5t4000) ( Rll t K) yKsKR,KE ) 

IFCIM .EQ. 1) REA0(5,400D) I WC 1 y K )tK=KR yKE ) 

IFCIU .EQ. 1) WRITE(6y2500) <U ( 1 y K) yK=l yKE > 

IFCIU .EQ. 1) HRITE(6y2500) I V 1 1 y K) yKsl yKE ) 

IFCIPT .EQ. 1) HRITE(6y2400) (PT ( lyK) yKsLyKE) 

IFIIRR .EQ. n WRITE(6y2500) (R U y K) y K^KRyKE ) 

IFCIU .EQ. 1) WRITE(6y2500) CMC 1 y K) ,K=1 yKE ) 

IFCIPT .EQ. 1» GO TO 45 
GO TO 60 

PT AND W FROM VARIABLE LIST 
45 CONTINUE 

INITIALIZE ONLY THE UPWIND COLUMN 

IFCIII .EQ. 2) GO TO 501 

INITIALIZATION ALL COLUMNS 

DO 50 I=^2yIE 
DO 50 K=1»KE 
UCl yK)-U(lyKl 
VCIvKIsVClyK) 

PT(IyK)=PTClyKl 

RCIyKIsRClyK) 

50 WCIyK)=HClfK) 

60 CONTINUE 

INITIALIZATION. EXCHANGE COEFFICIENT. INTEGRATED LIQUID WATER, AND 
SPECIFIC HEAT OF MOIST AIR 

501 CONTINUE 

CL=14.*G*.16/PTC1,1 ) 

lNTClyl)=0.0 

TClyll=PTCl,l) 

DO 80 K=2,KE 

PRCK)aC 1000. /PCK))*4. 286 
OZCK)sC 1.>ZAK/2.)«CDZACK)**2) 

TClyK)*PTCl.K)/PRCK» 

DW2=SQRTCCU(1,K»-UC1,K-1) )**2+CVCl.K»-VClyK-in**2)/DZACK) 

IFCHClyK) .GT. 0.0> GO TO 65 
CPTC1.K)=CP 

INITIALIZATION OF EXCHANGE COEFFICIENTS COVERING NEUTRAL, STABLE AND 
UNSTABLE CONDITIONS BETWEEN K=2 AND THE SURFACE 

S=C PTC lyK)-PT(l,K-l » I/DZACK) 

IFCK.EO. 2) SS=S 
GO TO 70 

65 CPTC l.K»=CP+CH*RSFCTClyKlyPCKn/CTCl,K)**21 
S*(TClyKJ-TClyK-l) »/DZACKI+CS/CPTCl,K) 

IFCK .EQ. 2) SSx^S 
70 IFCK .EQ. 2) GO TO 71 
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IFISS .GT. 0.) GO TO 75 

IF(SS .EQ. 0.) GO TO 77 

PHI=0SQRT(DSQRTU#D0-l5.DC^ZCKl/TLn 

KHfl«K)= .4*Z CK)*UF1*PHI 

KA(1,K)=KHC1»K)*PHI 

GO TO 79 

75 KNU«K)=( .4*Z(K)*UFl)/< 1 .+5 .♦ ( Z( K l/TL )) 

GO TO 78 

77 KMU,K)=.4*ZIK)*UF1 

78 KA(ltK)=KH(ltK) 

GO TO 79 

71 KM(X,2)=.16*$0RT(U( 1,2)*»2+VC lf2)**2 l♦DZA(2)/(AL0G(2A 

1(2)/Z0)»*2) 

KAU,2)=KM(l,2) 

UFl=SQRTCKM(l*2l*0HZ) 

IF(S .EO. 0.0) GO TO 8C 
TL=(UF1«3)*PT(1,2»/C.4*G*KA( 1,2 )*S) 

GO TO 80 

SHIR-TYPE EXPONENTIAL DECREASE 

79 DEC=eXPC-8.*F*ZUO/UFl) 

KM( I,K)=KH( 1,K»*DEC 

KA< 1,K)=KA( ltK)*DEC 

80 INT< 1,K) = INTU,K-1 H-CI*(W( 1,K)+^.667^W( 1 ,K-1 ) . 6671 *DZA (K ) 

CPT(1,1 )=CPT| 1,2) 

KA(l,l)=RA(l,2) 

IFCIII .EQ. 2) GO TO 502 
INITIALIZATION ALL COLUMNS 

00 81 1=2, IE 

DO 81 K=1,KE 

IF(K .EQ. 1) W(1,K)=0.0 

T(I,K)=TU,K) 

CPT(I,K)=CPT(1,K) 

KH(I,K)=KM(1,K) 

KAU,K)=KA(1,K) 

81 INTC I,K)=INT( 1,K) 

500 CONTINUE 

502 CONTINUE 

INITIALIZATION OF SURFACE TEMPERAIURE DIFFERENCE 

IFIIDTEM .EG. -1) GO TO 86 
DO 85 1=1, IN 

IFCIDTEM .EG. 0) DTEM(I)=0.0 
iFdOTEM .EG. 1) OTEMd )=DTEM1 

IF((IDTEM .EG. 1) .AND. (I .LT. ITEMD) OTEM(I)=O.C 
IFMIDTEH .EG. 1) .AND. (I .GT. ITEMR)) OTEH(I)=C.G 

85 CONTINUE 
GO TO 87 

86 REAO(5,1400) ( DTEM ( I ) , I ^2 , IN ) 

WRITE|6,2400) ( DTEM ( I ) , 1=2 , IN ) 

DTEM(1)=0.0 

DTEM(IE)=0»0 

UPPER BOUNDARY CONDITION FOR IMPLICIT INTEGRATION 

87 EU(KE)=0. 

EV(KE)*0« 
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EPT<KE):=0.0 

ER(KE)>0.0 

EW(KEI=0.0 

FUCKEI-UUfKE) 

FV(KE) = VU,KE) 

FPT(KE)=PT( 1,KE) 

FR(KE)aRilfKE) 

FWiKEjsHUf KE) 

OUTPUT TIME ,ENO TIME, AND TIME STEP CONTROL 

TIME=0,0 

PRT=OT 

90 IFCTIME .LT. PRT} GO TO 95 
PRT=TIME^0T 
GO TO 200 

95 IFITIME .GE. ET ) GO TO 10 
TIME=:^TIME+DT 
CALL STEP 
GO TO 90 
200 CONTINUE 


OUTPUT 


EDDY HEAT FLUX OUTPUT 

300 MRITE(6,3000) 

DO 311 1=1, IE 
DO 310 K=3,KN 

IF(CPTU,K) .LE. CP) HC{I,K)=-DEN»CP»(<KA(I,K^1)*(PTII,K+1)-PT(1,K 
in/(rA(K+ll-ZA(Kin*(KA(I,K)^(PT<I,K)-PT(I,K-lH/(ZA(K)^ZA(K-in I ) 
2*2. I 

310 IFCCPTII.K) .GT. CP) HC U ,K )=-DEN’K:PT ( I ,K) ♦ ( ( (KA ( I ,K) +KA ( I , K-H ) )/2 
!.)♦( <TCI,K+l)-TCIfK-l) )/UAU + l)-ZAlK-l ))+CS/CPT( I,K) ))*4.2 

WINDV2=S0RT(U| I ,2)**2^V(I,2)»*2)*.07 

IF(CPTCI,2) .LE. CP) HC(I,2)=-0EN*CP*60.*.16^WIN0V2*(PT(I,2)-PT(I, 
11) )/(ALQG(ZA(2)/ZQ)**2) 

IF(CPTII,2) .GT. CP) HCat2)=-DEN*CPT(I,2)*60.*.16*WlNDV2*(T(I,2)- 
lTa,l)^CS*ZA( 2)/CPT(I,2))/(AL0GUA(2)/Z0)»*2) 

HCCI,1)=HC(I,2) 

HCI1,KE>=HC(I,KN) 

311 CONTINUE 
CALL PRNT(HC,FMTH) 

C 

C LOCAL FRICTION VELOCITY OUTPUT 

C 

WRITEC6t3008) 

DO 313 I-1,IE 
DO 312 K=3,KN 

FLUX=(KM(I,K+1)»SQRTUU(I,K+1)-U(1,K) )*»2+ ( V( I ,K + 1 )-V (I ,K) ) ♦♦2 ) /DZ 
lA(K^l) ♦KM(I,K)*SQRTUU(1,K)-U(I,K-1))#»2 + (V(1,K)--V(I,K-1))*^2)/DZ 
2A(K))/2. 

IFIFLUX .GT. 0.) HCUfK)=SQRT(FLUX) 

312 CONTINUE 

WINDV2=SQRT(U( 1 ,2)**2+V< I ,2 ) ♦♦Z ) 

HCl I,2)=.A*WINDV2/AL0G( ZA(2)/Z0) 

HC(l,l)=HC(I,2) 

313 HCCI,KE)=HC(I,KN) 
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CALL PRNTIHCf FMTH) 


C 

C 

C PRINT AND PUNCH OF PROGNOSTIC VARIABLES FOR INPUT AS UPWIND CONDITIONS 

C 

C 

C WIND OUTPUT 

C 

WRITEI6«3009) 

DO 314 I-ltlE 
00 314 K«ltKE 

314 HC(lrK)sl.0E-2*U(ltKl 
CALL PRNTIHCfFMTHJ 
WRITE(6t30lC) 

DO 316 I=1>1E 
00 316 K=1*KE 
316 HCI1.K)=1.0E-2*V(I.K) 

CALL PRNTIHCtFMTH) 

WRITE(7t4000) ( U( IZ,K ) , K=1 , KE ) 

WRITE(7,4000) ( V( IZ«K) >K^ltKE) 

C 

C TEMPERATURE OUTPUT 

C 

WRITE(6»3100) 

CALL PRNT(TfFMTT) 

WRITEI7,5000) ( PTC IZ ,K > *K=1 .KE ) 

C 

C MIXING RATIO OUTPUT 

C 

WRITEC6t3200) 

CALL PRNTCR,FMTH) 

WRITEC7.4000) ( R ( IZ * K) t K=1 ,KE ) 

C 

C LIQUID WATER CONTENT OUTPUT 

C 

WRITE! 6.3300) 

DO 315 1=1, IE 

00 315 K=1,KE 

315 HC(I,K)=DEN*WCI,Kl*1.0E+6 
CALL PRNTCHC.FMTH) 

WRITE(7,4000) ( W( IZ.K) , K=1 ,KE ) 

C 

C DEW POINT DEPRESSION OUTPUT 

C 

WRITE(6,3400) 

DO 320 1 = 1. IE 
DO 320 K=1.KE 
C 

C DEW POINT DEPRESSION COMPUTATION 

C 

E=P(K)*R(I.K)/( .62465«R(I.K) I 
Tl=T(I,K) 

DO 316 M=I,3 
RS=RSFCT1,P(K)) 

ES=P(K)*RS/« .62465+RS) 

EE=E/ES 

318 Tl»Tl*(Tl-35.a6)4AL0GCEE)/17.26939 
320 HCII.K)=T(I,K)-T1 
CALL PRNTCHC.FMTH) 

C 

C TURBULENT EXCHANGE COEFFICIENT OUTPUT 
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WRITE(6«350C) 

DO 331 1=1, IE 
DO 330 K=3,KN 
EX=IKA( 1,K)^KA( l,K*l))/2. 

330 HC(I,K)=EX 

HC<I,2)=.16*$QRT<U( I,2)**2+V(I,2 )**2 ) *2 A( 2 ) /ALOG ( Z A( 2 )/Z0) 

HCi I,1)=G. 

331 HC(1,KE)=HC(1,KN) 

CALL PRNT(HC,FMTH) 

IFCIP .EQ. 0) GO TO 95 

RADIATIVE FLUX OUTPUT 

WRITE(6,3600) 

DO 3^ 1=1, IE 
DO 340 K=1,KE 

HR=RF*S1GMA»(T( I,1)**4)*EXP(“INTI I,KE)+INT( I,K) )*60. 

340 HCI1,K)=HR 

CALL PRNT(HC,FMTHJ 

RADIATIVE COOLING RATE OUTPUT 

WRITE(6,3700J 
DO 350 1=1, IE 
DO 350 K=1,KE 

C00L=-3600.*( Wil,Kl^*.667)*CR*(T( I , 1 ) **4. ) *EXP ( -I NT ( I ,KE ) + INT ( 1 , K ) 
1)*CP/CPT(I,K) 

350 HCU,K)=COOL 

CALL PRNT(HC,FMTH) 

WATER VAPOR FLUX OUTPUT 

WRITE(6,3800) 

DO 361 1=1, IE 
DO 360 K=3,KN 

HCCI,K)=-DEN*.5*(<KA(I,K+1)*(R( I,K+1)-R(I,K))/(ZA(K^1)-ZA<K)))+(KA 
1(I,K)*(R(I,K)'RCI,K-1) )/(ZA(K)-ZA(K-D) )) 

360 CONTINUE 

HINDV2=SQRT(UU,2)**2+V(I,2)**2) 

HCI I,2)=-DEN*.lt>*WINDV2*IRiI,2)-R<I , 1 ) ) /( ALOG I ZA ( 2 )/20) ♦♦2 ) 
HC(I,1)=HC(I,2) 

HC(I,K£)=HC( 1,KN) 

361 CONTINUE 

CALL PRNT(HC,FMTH) 

GO TO 95 
400 STOP 

1000 FORMAT! 10E8 .2) 

1100 FORMAT! 1415) 

1200 FORMAT(3E10.3,615,-2PF6.0) 

1205 FORMAT! 11A4/11A-,) 

1300 FQRMAT(4E10.3) 

1400 F0RMAT(8F10.2) 

1500 FORMAT! 8E10. 3) 

2000 FORMAT! 1H1,////,60X,10HINPUT DATA,//, IHO, IPIOEIO .3 ) 

2100 FORMAT! IHO, 1415) 

2200 FORMAT! IHO, 1P3E 10. 3 , QP6 15 ,-2PF6 . 0 ) 

2300 FORMAT! IHO, 4E1C. 3) 

2400 FORMAT! IHO, 8F10.2) 

2500 FORMAT! IHO, 6E10. 3) 
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3000 FORMAT ( 30H1 EDDY HEAT FLUX IM WATTS/ CM** 2/// I 
300B F0RMAT(28H1FR1CTI0N VELOCITY IN CM/SEC///) 

3009 FORMAT! 29H1U COMPONENT OF WIND IN M/SEC///I 

3010 F0RMAT(29H1V COMPONENT OF WIND IN M/SEC///I 
3X00 FORMAT! 21HXTEMPERATURE IN DEG K///J 

3200 FORMAT! 20HIMIXING RATIO IN G/G///I 

3300 FORMAT! 31H1LI0UID WATER CONTENT IN G/M**3///t 

3900 FORMAT !3CHIDEH POINT DEPRESSION IN DEG C///J 

3500 FORMAT !49M1 TURBULENT EXCHANGE COEFFICIENT IN CM**2/SEC///) 

3600 FORMAT! 31H1RAD1ATIVE FLUX IN CAL/CM**2M IN/// » 

3700 FORMAT! 30H1RAD1ATIVE COOLING IN DEG C/HR/// ) 

3800 FORMAT !39HI EDDY MOISTURE FLUX IN G/CM**2 SEC///) 

C 

C 2-FORHAT EQUALS HEXADECIMAL FORMAT 

C 

4000 FORMAT !10Z8) 

5000 FORMAT! 5216) 

END 

C 

C INTEGRATION SUBROUTINE 

C 

C INTEGRATE ONE TIME STEP AND COMPUTE NEW PROGNOSTIC C DIAGNOSTIC 

C VARIABLES 

C 

SUBROUTINE STEP 

REAL*8 PT,EPT,FPT,DPT, AT,CT,BT 
REAL*4 KAfINTtKW«LtLM,LMAX,KM 

COMMON PT!40*6G ),EPT!60) ,FPT!60) »DPT«AT .CTtOTt 
a R140>6u) tW(40«60)tU!40«60) *T!40f 60) *KA!40t60) * 

7KM1 40t60) , 

1 lNT(40»60),CPTl40t60) «HC!4C,60)tX!40)«OX!40),DTEM!40), 

2 P(60)tZA!60)tOZA!60)t ER! 60 ) tFR 160) »EW! 60 ) t 

3 FW!6C) tEU!60) tFU! 60 ) t Z!60 ) .DZ! 60 ) ,PR !60 ) * L f DEN tCPtGt RAtRWt 

4 SIGMA,TIME»DT,TIH,UF,RF»KWtZOtCVfUItDTEMI,ZAK,XAItDELXfCC« 

5 CH»CIfCK.CL,CR.CS>UU>FV!60),EV!60)«V!40t60)t F» 

6 KEfKNtlEf lN>lL«lR,ISED.IRAO,lRSFCtIOTGHtITEMLiITEMR, IP, II 
C 

C UPWIND DO LOOP OVER HORIZONTAL GRID SYSTEM 

C 

DO 80 M=2,IN 
I=IN*2-M 
C 

C DOWNWARD DO LOOP OVER VERTICAL GRID SYSTEM TO SET UP IMPLICIT 

C INTEGRATION 

C 

00 20 N=2yKN 

K*KN+2-N 

OD=DT/DZ!K) 

A=DD*KA!I,K*1 )/!X.-«^ZAK) 

AWsA 

AU*D0*XM!I,K*1)/!1.+ZAK) 

AUU-II*AU 

AT=A 

C3D0*KA!I,K) 

CW=C 

CU3D0*KM!I,K) 

CUU=II*CU 

CT=C 

IF!ISED .NE.l) GO TO 10 
C 

C DROP SEDIMENTATION 
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I III 


c 

0SsCV*0D*0ZA(K)/2. 

AW«AW^0S«(W(ItK«l)*«.667} 

CWsCW-0S*(W«ItK-l)**.667l 

10 B*l.*A+C 

BUU-1.-»AUU-*’CUU 
BT=l.*AT*CT 

DPT=PTIItK)-U<I»K)*(P7(I»K)-PT(I-l,K))*OT/OX(I» 

RAOIATIOMAL COOLING 

IFUIRAO .EQ. 1) .AND. «WII,K» ,GT. 0.0)) DPTaOPT-CR*Wt I,K )*(PT Cl. 
n)**4.)*EXP(-INTII,KE)*INTCI»K))*PR(K)*0T/CHCI,K)**.333) 

0U=U(I»K)*(1. + (II-1)*(AU*CU) )-( II“l)*IAU*UII»K+l XU*U(l»K-l) ) 
1-OCI»K)*IU<I,K)-UCI-1,K))*DT/OXCI) ♦ F*IVC I ,K)-V( I »KE ) )*0T 

DV=VCI.K)*Cl.+CII-l)*CAU+CU))-< II-D* C AU*VC I.K^l )+CU*VC I.K-1 ) ) 
1-U(I,K)*CV(I,K)-VCI-1,K ))*DT/DXCI) ♦ F*IU( I.KE )-UC I »K) )*DT 
OH = W<I.K)-UCI«K)»C WCI.K)- W(I-1.K))*0T/0X(I) 

OR = R(I.K)-UC I.K)*C RCI.K)- R( I-l.K) )*0T/DXI I ) 

EUC K )sCUU/C BUU-AUU*EUC K^l ) ) 

EVi K)-CUU/C BUU-AUU*EVC K+l ) ) 

EPTCK)=CT/CBT-AT*EPT(K+1) ) 

ER(K)«C/IB-A*ERCK41)) 

ew(k)>ch/cb-ah*ewck4'1) ) 

FU(K)s(OU'»AUU*FU(K'»l))*EUCK)/CUU 
FVCK)=IDV+AUU*FV(K+1))*EVCK)/CUU 
FPT(K)=(OPT*AT*FPT(K+l) )*EPT(K)/CT 
FRlK)«(DR+A*FRIK*l))*ERCK)/C 
FW(K)*(0W-»AW*FHCK*1) )*EMIK)/CW 
20 CONTINUE 

UPDATE SURFACE BOUNDARY CONDITION 

IFCTINE .LE. TIM) PTC I . l)sPT ( 1. 1 ) +0TEMC I )*TIME/TIM 
IFCTIME .LE. TIM) T( I , I )=PT C 1 . 1) 

IFIIRSFC .EQ. 0) RII.l)sFRC2)/(l.-ER(2) ) 

SATURATION SURFACE BOUNDARY CONDITION ON R 

IFCIRSFC .EQ. I) RCI.D^RSFITCI.D.IOOO.) 

UPWARD DO LOOP OVER VERTICAL GRID SYSTEM TO COMPUTE NEW PROGNOSTIC 
AND DIAGNOSTIC VARIABLES 

INTI I.1)=0.0 
DO 70 K=2,KE 
IFCK .EQ. KE) GO TO AO 
UCI.K)*EUIK)*U4I,K-1)+FUCK) 

V(I,K)»EVCK)*VCI.K-1)+FVCK) 

PT(I.K)=EPTCK)*PTCI»K-1 )*FPTIK) 

R(I,K)»ERCK)*R( IfK-D + FRCK) 

H(I,K)=EWIK)*W(I,K-l)+FH(K) 

COMPUTE TEMPERATURE 

TCI .K)>PT(I.K)/PR(K) 

SATURATION ADJUSTMENT 
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RSsRSF(T(ItK>*PIKn 

IP((R(ItK) .LE. RS) .AND. (WCItK) .LE. C.O) ) GO TO 40 
DR»(RIIfK»-RS)/<l.+CC*RS/«T( I»K»**2J) 

IF(RdtK) .GT. RSI GO TO 30 
IF(-0R .LE. W(ItKI) GO TO 30 
DR=-W(1 tKI 

30 TdtKl»Tn.K)+DR*L/CP 
PT(ItK)sT(I»KI*PR(K) 

R(ItK}«R(I«K>-OR 

H(I»K)-M(ItK)4DR 

COMPUTE NEW VALUES OF INT. CPT. AND KA 

40 IFIWdtK) .LT. 0.0) W(I.K)=0.0 

INTUtKI^INTU tK-1 >'«-CI*(Wd»KI««.6674W( ltK-l)**.667)*DZA(K I 
DW2*SQRTnUd»K)-Ud,K-in*»2+<Vd.K)-VdtR-l))**2)/DZA(K) 
IF(Wd.K) .GT. 0.0) GO TO 50 
CPT)I»K)*CP 

DIAGNOSIS OF EXCHANGE COEFFICIENTS COVERING NEUTRAL. STABLE AND 
UNSTABLE CONDITIONS BETWEEN K^2 AND THE SURFACE 

S*<PTd.K)-PT( I,K-1) )/DZA(K) 

IFIK.EO. 2) SS=S 
GO TO 60 

50 CPTCI.K )=CP+CH*RSF(Td ,KI,P(K) )/(T( I,K)**2) 

S=iTd.K)-7 d,<-l) )/DZACK)+CS/CPT(I,K) 

IFIK .EO. 2) SS=S 
60 IF(K .Eu. 2) GO TO 71 
IF(SS .GT. 0.) OO TO 75 
IF(SS .EQ. 0.) GO TO 77 
PHI =DSQRT(DSQRTd. DO-15. D0*Z(K)/TL) ) 

KM(I.K)s .4*Z (K)»UF1*PHI 
KA( I,K)sKM( ItK)*PHI 
GO TO 79 

75 KH( I.K)s< .442«K>*UF1)/(1.45.*»ZIK)/TL)) 

GO TO 78 

77 KM(I»K)s.4*Z(K)*UFl 

78 KA( I«K)=KM( ItK) 

GO TO 79 

71 KMdt2)=.16*SQRT(U( 1 .2 )>**2 4V( I ,2 )«*2 ) *OZA ( 2 )/ ( ALOG ( ZA 

H2)/Z0)**2) 

KA( l«2)sKHdt2) 

UFl=SQRT(KHd t2)*DWZ) 

IFIS .EQ. 0.0) GO TO 70 
TL»)UFl«3)*PTd.2)/«.4*6*KA(I,2 )*S) 

GO TO 70 

SHIR-TYPE exponential DECREASE 

79 DEC*EXP(-8.*F*Z(K)/UF1) 

KM(IiK)-KM(ItK)«DEC 

KA( ItK)sKA( I.K)*DEC 

70 CONTINUE 

UPDATE SURFACE BOUNDARY CONDITION ON R AFTER SATURATION ADJUSTMENT 
IF NO FLUX CONDITION USED 

IFdRSFC .EQ. 0) R(I.l|BRdt2) 
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80 CONTINUE 
C 

C UPDATE DOWNWIND BOUNDARY CONDITION 

C 

DO 90 K=X,KE 
PT(IEvK)sPT(INtK) 

U|IE»K)>UIIN»K) 

V(IE,K)=VCIN»K) 

R(IEtK)=RilNtKl 

W(IEtK)=W(INtK| 

KA( IE«K l-KAdN.K) 

KH( lEtK)s:KM(INtK) 

T4IEtK)=T(IN,K) 

CPT|IEtK)=CPT(IN,K) 
iNTCIEf K) = INTaN,K) 

90 CONTINUE 
RETURN 
END 
C 
C 
C 

C PRINT SUBROUTINE 

C 

SUBROUTINE PRNT(OfFORM) 

REALMS PT*EPT»FPT,DPTfAT,CT.BT 
REAL«4 KA»INTtKW»L»LM»LMAX»KM 

COMMON PT(40f 60),EPT(60)tEPT(60) t DPT, AT *CT , BT, 

8 R(40,60),W(40,60),U(40,60) ,T(40,60) ,KA (40,601 , 

7KM( 40,60), 

1 INT(40,60),CPT(40,60I,HC(4C,60),X(40|,DX(40),DTEM(40), 

Z P(60),ZA(60),0ZA(6C), ER( 60 1 ,FR (601 , EW( 60 ) , 

3 FWC60) ,EU(60) ,FU(60),Z(60)tDZ(60),PR(6C),L,DENrCP,G,RA,RW, 

4 SIGMA, TIME, OT,TIM,UF,RF,KW,ZO,CV,UI,DTEMI,ZAK,XAI,DELX,CC, 

3 CH,CI,CK,CL,CR,CS,UU,FV(6C),EV(60),V(40,60), F, 

6 KE,KN,IE,1N,IL,IR,ISED,IRAD,IRSFC,IDTEH,ITEML,1TEHR,1P,II 

DIMENSION 0(40, 601, FORM(ll) 

C 

C OUTPUT DOCUMENTATION 

C 

WRITE(6,4000) T IME ,OTEM I , ITEML, I TEMR 

C 

C PAGE 1 COLUMNS 1-10 

C 

WRITE(6,4100) (X(1),I«I,10> 

00 430 J=1,KE 
K=KE-J+1 

450 WRITE(6,F0RM) ZA(K ) , (0 ( 1 ,K) , 1=1, 10) 

C 

IF( IE .LT, 11) GO TO 490 
C 

C PAGE 2 COLUMNS 11-20 

C 

WRITE(6,4103) ( X( I ) , 1=1 1 ,20 ) 

00 460 J=1,KE 
K=KE-J^1 

460 WRITE(6,F0RM) ZA(K ) , (0 ( I ,K ) , 1=11 , 20 ) 

C 

IF(IE •LT. 21) GO TO 49(J 
C 

C PAGE 3 COLUMNS 21-30 

C 
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WRlTE(6t4l05) (XU) tl=21»30) 

DO 470 J=lfKE 
K«KE-J+1 

470 WRITE( 6 »F 0 RM) Z A(K ) » (0 ( I ^K) t 1^21 « 30 ) 

C 

IFCIE *LT« 3U GO TO 490 

PAGE 4 COLUMNS 31-1 E 

WRITE(6«4X03) ( X( I ) t Is31« IE ) 

DO 480 J^ItKE 
K=KE-J+1 

480 WRlTE( 6 fF 0 RH) Z A(K) * (0( 1 tK) t I =31 t IE ) 

490 CONTINUE 

4000 FORMAT(6H0TlHE=fF7.0t4H SEC 1 6 X 96 HOT EMI = tF5 • 1 1 2M K » 6 X » 6 HITEML= 1 1 2 12 
IXt6HITEMR=*I2) 

4100 FORMAT! lH0«6HX(KH)st6X9 10(-^PFI2.2) / IH «3Xf5HZ(CM)) 

4103 FORMAT! IHl«6HX(KM)=,6X»I0(-3PF12.2) / IH »3Xt5HZ(CM)) 

RETURN 

END 


SATURATION MIXING RATIO AS A FUNCTION OF TEMPERATURE A AND 
PRESSURE B 

FUNCTION RSF(AvB) 

ES* 6 . 1078*EXP (17.26939* ( A-273. 16 )/(A-35 , 86 ) ) 
RSF=.62465*ES/CB-ES) 

RETURN 

END 


DATA CARDS 

• 20E^2 I.80E^3 3.6CE^3 1.80E4>3 

2.50E-1 

•71E^1 

2 

.87E+2 l.OOE-2 O.OOE+0 

5.C0E-»-C 

I I 1 1 1 1 

1 1 

1 

1 

1 1 43 9 


l.OOOE^l 0 . 200 E -«^0 2.000E4^0 

(IH 1PE11.2910E12.2) 

(LH 1PE11.2«0PI0F12.2) 

10 43 

1 

6 

2 9 1000. 

FMTH 

FMTT 

1.090E-2 2.870E+2 

PROGNOSTIC VARI ABLE(UtV«PT 

«R AND W) 

INPUT 

IF 

USED 
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